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'^^ ' Abstract 

We derive some new structural results for the transfer matrix of square- 
. lattice Potts models with free and cylindrical boundary conditions. In partic- 

f— ^ \ ular, we obtain explicit closed-form expressions for the dominant (at large |g|) 

diagonal entry in the transfer matrix, for arbitrary widths m, as the solution 
- of a special one-dimensional polymer model. We also obtain the large-g expan- 

r> \ sion of the bulk and surface (resp. corner) free energies for the zero-temperature 

C^ ' antiferromagnet (= chromatic polynomial) through order g~^^ (resp. q~^^). Fi- 

nally, we compute chromatic roots for strips of widths 9 < m < 12 with free 
boundary conditions and locate roughly the limiting curves. 
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1 Introduction 

The Potts model [1-4] on a regular lattice C is characterized by two parameters: 
the number q of Potts spin states, and the nearest-neighbor coupling v = e^'^ — 1} 
Initially g is a positive integer and f is a real number in the interval [—1, +cxd), but 
the Fortuin-Kasteleyn representation (reviewed in Section 2.1 below) shows that the 
partition function Zciq, v) of the q-staXe Potts model on any finite graph G is in fact 
a polynomial in q and v. This allows us to interpret q and v as taking arbitrary real 
or even complex values, and to study the phase diagram of the Potts model in the 
real (g, f )-plane or in complex (g, t>)-space. 

According to the Yang-Lee picture of phase transitions [5], information about 
the possible loci of phase transitions can be obtained by investigating the zeros of 
the partition function for finite subsets of the lattice C when one or more physical 
parameters (e.g. temperature or magnetic field) are allowed to take complex values; 
the accumulation points of these zeros in the infinite-volume limit constitute the 
(possible) phase boundaries. For the Potts model, therefore, by studying the zeros of 
Zc^q^v) in complex (g,f)-space for larger and larger pieces of the lattice £, we can 
learn about the phase diagram of the Potts model in the real (g, w)-plane and more 
generally in complex (g,f)-space. 

The partition function for mxn lattices can be efficiently computed using transfer 
matrices. Though the dimension of the transfer matrix (and thus the computational 
complexity) grows exponentially in the width m — thereby restricting us in practice 
to widths m < 10-30 — it is straightforward, by iterating the transfer matrix, to 
handle quite large lengths n. Indeed, by implementing the transfer-matrix method 
symbolically (i.e., as polynomials in g and/or v) and using the Beraha-Kahane- Weiss 
theorem (reviewed in Section 2.3), we can handle directly the limit n — *> oo and 
compute the limiting curves Bm of partition-function zeros. At a second stage we 
attempt to extrapolate these curves to m = oo. 

Since the problem of computing the phase diagram in complex (g, f )-space is dif- 
ficult, it has proven convenient to study first certain "slices" through (g, f)-space, in 
which one parameter is fixed (usually at a real value) while the remaining parameter 
is allowed to vary in the complex plane. One very interesting special case is the chro- 
matic polynomial {v = —1), which corresponds to the zero-temperature limit of the 
Potts antiferromagnet (/3J = — oo). In previous papers [6-10] we have used symbolic 
transfer-matrix methods to study the square-lattice and triangular-lattice chromatic 
polynomials for free, cylindrical, cyclic and toroidal boundary conditions.^ '^ In this 



"'^Here we are considering only the isotropic model, in which each nearest-neighbor edge is assigned 
the same coupling v. In a more refined analysis, one could put (for example) different couplings 
wi,W2 on the horizontal and vertical edges of the square lattice, different couplings wi,W2,i'3 on the 
three orientations of edges of the triangular or hexagonal lattice, etc. 

^See also the bibliographies of [6-10] for reference to the important related works of Shrock and 
collaborators. 

■^We adopt Shrock's [11] terminology for boundary conditions: free (m-p x np), cylindrical (mp x 
np), cyclic (mp x rip), toroidal (mp xnp), Mobius (rrip x n^p) and Klein bottle {mp x nxp)- Here the 



paper we provide some new structural results for the transfer matrices with free and 
cylindrical boundary conditions. (For simplicity we restrict attention to the square 
lattice, but the methods could easily be adapted to handle the triangular lattice.) 
In particular, we shall obtain explicit closed-form expressions for the dominant (at 
large |g|) diagonal entry in the transfer matrix, for arbitrary widths m, by solving a 
special one-dimensional polymer gas. We shall also obtain similar but weaker results 
for the dominant eigenvalue. Finally, we shall obtain the large-g expansion of the 
bulk free energy through order g~^^, extending by 11 terms the previous computation 
of Bakaev and Kabanovich [12], and the large-g expansions of the surface and corner 
free energies through orders g~^^ and g"^^, respectively. 

This paper is organized as follows: In Section 2 we review the Fortuin-Kasteleyn 
representation [13, 14] of the g-state Potts model, the basic facts about transfer ma- 
trices in this representation [6,15], and the Beraha-Kahane-Weiss theorem [16-20]. 
In Section 3 we prove some new structural properties of the transfer matrix for a 
square-lattice strip of width m and free or cylindrical boundary conditions; in partic- 
ular, we obtain closed-form expressions for the dominant entry of the transfer matrix, 
for arbitrary widths m. In Section 4 we study the large-g expansion of the leading 
eigenvalue of the transfer matrix; our results are similar to those obtained for the 
dominant entry, but less explicit. In Section 5 we study the limit m -^ oo of the 
strip free energy; among other things, we compute the large-g expansions of the bulk, 
surface and corner free energies, and we carry out a differential-approximant analysis 
to locate the singularities of those free energies in the complex g-plane. In Section 6 
we provide some additional information concerning the chromatic roots of strips of 
widths 9 < m < 12 with free boundary conditions. Finally, in Section 7 we list some 
open problems for future research. In Appendix A we prove some identities arising 
in the study of the dominant transfer-matrix entry for both boundary conditions. In 
Appendix B we discuss a conjecture concerning the upper zero-free interval for real 
chromatic roots of bipartite planar graphs. 

2 Preliminaries 

In this section we review briefly some needed background on chromatic and Tutte 
polynomials (Section 2.1), transfer matrices (Section 2.2) and the Beraha-Kahane- 
Weiss theorem (Section 2.3). 

2.1 Chromatic polynomials, Potts models, and all that 

Let G = {V, E) be a finite undirected graph, and let g be a positive integer. Then 
the qr-state Potts-model partition function for the graph G is defined by the 



first dimension [m) corresponds to the transverse ("short") direction, while the second dimension 
(n) corresponds to the longitudinal ("long") direction. The subscripts F, P and TP denote free, 
periodic and twisted-periodic boundary conditions, respectively. 



Hamiltonian 

^Potts(cr) = - ^ Je5{ai,aj) , (2.1) 

e=ijGE 

where the spins cr = {(Xjjjgy take values in {1,2,..., q}, the Jg are couphng constants, 
and the 6 is the Kronecker delta 

r/ ,x f 1 if a = 6 , ^, 

The partition function can then be written as 

^S°"^(g,v) = Yl n [l + ^^'^(^-^^)]' (2-3) 

0-: y-+{l,2,...,g} e=ijeE 

where Ve = e^'^^ — 1. Please note, in particular, that if we set Ve = —1 for all 
edges e, then Zq°^^^ gives weight 1 to each proper coloring and weight to each 
improper coloring, and so counts the proper colorings. Proper g-colorings (fg = —1) 
thus correspond to the zero-temperature {f3 — * +oo) limit of the antiferromagnetic 
(Je < 0) Potts model. 

It is far from obvious that Z^°**^(g, v), which is defined separately for each positive 
integer q, is in fact the restriction to g G Z_|_ of a polynomial in q. But this is in fact 
the case, and indeed we have: 

Theorem 2.1 (Fortuin Kasteleyn [13,14] representation of the Potts model) 

For integer q > I, 

ACE eGA 

where k{A) denotes the number of connected components in the subgraph {V,A). 

Proof. In (2.3), expand out the product over e E E, and let A C E be the set 
of edges for which the term VeS{ai,aj) is taken. Now perform the sum over maps 
a: V ^ {1,2, . . . ,q}: in each component of the subgraph {V, A) the color (Tj must be 
constant, and there are no other constraints. We immediately obtain (2.4). D 

Historical Remark. The subgraph expansion (2.4) was discovered by Birkhoff 
[21] and Whitney [22] for the special case Ve = —1 (see also Tutte [23,24]); in its 
general form it is due to Fortuin and Kasteleyn [13,14] (see also [25]). 

The foregoing considerations motivate defining the multivariate Tutte polyno- 
mial of the graph G: 

^G(g,v) = J^g'^^^^n^- (2.5) 

ACE e&A 

where q and v = {ve\e<^E are commuting indeterminates. If we set all the edge weights 
Ve equal to the same value v, we obtain a two- variable polynomial that is equivalent to 



the standard Tutte polynomial TG{x,y) after a simple change of variables (see [26]). 
If we set all the edge weights Ve equal to —1, we obtain the chromatic polynomial 

PG{q) = ZGiq,-l). 

Further information on the multivariate Tutte polynomial Zclq, v) can be found 
in a recent survey article [26] . 

2.2 Transfer matrices 

For any family of graphs Gn = {Vn, En) consisting of n identical "layers" with 
identical connections between adjacent layers, the multivariate Tutte polynomials of 
the Gn (with edge weights likewise repeated from layer to layer) can be written in 
terms of a transfer matrix [6, 15]. Here we briefly summarize the needed formalism [6] 
specialized to the case of an m x n square lattice with free and cylindrical boundary 
conditions. As usual, free (resp. cylindrical) boundary conditions means free (resp. 
periodic) boundary conditions in the transverse ("short") direction and free boundary 
conditions in the longitudinal ( "long" ) direction. These lattices are denoted by mp x 
Up and mp x np, respectively. 

Consider the m x n square grid with edge weights fi,i+i on the horizontal edges 
{1 < i < m) and Vi on the vertical edges {1 < i < m). Site m + 1 is always to 
be understood as a synonym for site 1. If the weight Vm,m+i = Vm,i is zero (resp. 
nonzero) we are considering free (resp. periodic) transverse boundary conditions. 

We fix the "width" m and consider the family of graphs G„ obtained by varying 
the "length" n; our goal is to calculate the multivariate Tutte polynomials ^^^(g, v) 
for this family by building up the graph Gn layer by layer. What makes this a bit 
tricky is the nonlocality of the factor g'^(^) in (2.5). At the end we will need to know 
the number of connected components in the subgraph (Vn, A); in order to be able to 
compute this, we shall keep track, as we go along, of which sites in the current "top" 
layer are connected to which other sites in that layer by a path of occupied edges 
(i.e. edges of A) in lower layers. Thus, we shall work in the basis of connectivities of 
the top layer, whose basis elements e-p are indexed by partitions V of the single-layer 
vertex set {1, . . . , m}. The elementary operators we shall need are: 

• The join operators 

■iijG-p = e-p,ij , (2.6) 

where V • ij is the partition obtained from V by amalgamating the blocks 
containing i and j (if they were not already in the same block). Note that all 
these operators commute. 

• The detach operators 

D.ep = <^ (2.7) 

[qep ii{i}eV 

where V\i is the partition obtained from V by detaching i from its block (and 
thus making it a singleton). Note that these operators commute as well. 

6 



Note, finally, that D^ commutes with Jjj whenever k ^ {hj}- 

The horizontal transfer matrix, which adds a row of horizontal edges, depends on 
the boundary conditions, and is given by 



m— 1 



Hf = n (1 + ^^.^+i-'m+i) (2-8a) 



i=l 



Hp = (l + t;„.,iJ„.,i)HF (2.8b) 

(note that all the operators in both products commute). The vertical transfer matrix, 
which adds a new row of sites along with the corresponding vertical edges, is 

m 

V = l[{vJ + D,) (2.9) 

(note once again that all the operators commute). The multivariate Tutte polynomial 
for Gn is then given [6] by the formula 

ZG„(g,v) = a.TH(VH)"-ieid , (2.10) 

where "id" denotes the partition in which each site i G {1, . . . , m} is a singleton, and 
the "end vector" oj^ is defined by 

uj'^ep = gl^l . (2.11) 

The transfer matrix is thus 

T = VH . (2.12) 

In principle we are working here in the space spanned by the basis vectors e-p 
for all partitions V of {1, . . . ,m}; the dimension of this space is given by the Bell 
number Bm [27-30]. However, it is easy to see, on topological grounds (thanks to 
the planarity of the Gn), that only non-crossing partitions can arise. (A partition is 
said to be non-crossing iia<b<c<d with a, c in the same block and b, d in the 
same block imply that a, 6, c, d are all in the same block.) The number of non-crossing 
partitions of {1, ... , m} is given by the Catalan number [28, 30] 

C.„ = ';"'" „ = ^- f2'») . (2.13) 

When the horizontal couplings fi,i+i are all equal to or —1 (which is the case for 
the chromatic polynomial with either free or periodic transverse boundary conditions), 
then the horizontal operator H is a projection (i.e., H^ = H), so that our vector space 
V splits up as a direct sum V = Vo © Vi, where Ht> = for t> G Vo and Hv = v for 
f G Vi. Then every vector f G Vq is an eigenvector of T = VH with eigenvalue 0; and 
for each eigenpair (A,f) of T = VH with v ^ Vq, the pair (A, Hv) is an eigenpair of 
T' = HVH. In this situation, therefore, we can work in the smaller vector space Vi 
by using the modified transfer matrix T' = HVH in place of T = VH, and using the 
basis vectors 

fp = Hep (2.14) 



in place of e-p. Indeed, if Te-p = '^^, t-ppre-p/ , then T'fp = ^p, t-p-p/fp', as is easily 
seen by applying H to both sides and using H^ = H. Please note that fp = if P has 
any pair of nearest neighbors in the same block. We thus work in the space spanned 
by the basis vectors fp where P is a non-crossing non-nearest-neighbor partition of 
{l,...,'m}. The dimensionality of this space depends on the transverse boundary 
conditions: 

• Free transverse boundary conditions: The number of non-crossing non-nearest- 
neighbor partitions of {l,...,m} is given [31, Proposition 3.6] [32] by the 
Motzkin number M^-i, where [28,30,33-36]*^ 

Ln/2J 

• Periodic transverse boundary conditions: The number of non-crossing non- 
nearest-neighbor partitions of {1, ... , m} when it is considered periodically (i.e. 
when 1 and m also are considered to be nearest neighbors) is given by [36, Sec- 
tion 3.2, family R2] 

I 1 for m = 1 , , 

dm = l^ f ^9 (2.16) 

I Km tor m > 2 

where the Riordan numbers (or Motzkin alternating sums) Rm [33,34,36]^ are 
defined by i?o = 1, -Ri = and 

m— 1 

R^ = ^(-l)'"-'=-iMfc form>2 (2.17) 

A:=0 

Finally, spatial symmetries further restrict the subspace whenever the couplings 
fj^i+i and Vi are invariant under the symmetry. Again the symmetries depend on the 
transverse boundary conditions: 

• Free transverse boundary conditions: Here the relevant symmetry is reflection 
in the center line of the strip. For reflection-invariant couplings, we can work in 
the space of equivalence classes of non-crossing non-nearest-neighbor partitions 
modulo reflection. The dimension SqFree(m) of the transfer matrix is then given 
by the number of these equivalence classes. The exact expression for SqFree(m) 
was obtained in [37, Theorem 2]: 

SqF.ee(™) = -M,,.-. + <^ g JJ,^^^;^, (2.18) 



^ Warning: Several references use the notation m„ to denote what we call M„ ; and one reference 
[34] writes M„ to denote a different sequence. 

'"'In most of the literature (e.g. [33,34]) these numbers are called 7^. We have adopted the recent 
proposal of Bcrnhart [36] to name them after Riordan [33] and denote them Rm ■ 



where 



m + 1 



m 



(2.19) 



and |_pj stands for the largest integer < p.^ The asymptotic behavior of 
SqFree(m) is given by [37, Corollary 1] 



SqFree(m) 



4v^ 



3^^-3/2 



l+0(- 

m 



as m ^ oo 



(2.20) 



Periodic transverse boundary conditions: For the square lattice with periodic 
transverse boundary conditions, both reflections and translations are symme- 
tries. We therefore define equivalence classes of non-crossing non- nearest-neighbor 
partitions modulo reflections and translations and the corresponding number 
SqCyl(m) of equivalence classes. To our knowledge there is no known closed 
form for these numbers; but there is a conjectured formula [39, Conjecture 2.2] 



\ [TriCyl(m) + iA^pp (f )] for even m 

SqCyl(m) 

for odd m > 3 

(2.21) 
where NY-p{m) is the number of different eigenvalues for a strip of either square 
or triangular lattice with cyclic boundary conditions (i.e., free transverse and pe- 
riodic longitudinal boundary conditions), and TriCyl(m) is the number of equiv- 
alence classes of non-crossing non-nearest-neighbor partitions modulo transla- 
tions. It is known [40] that 



iVFp(m) 



2(m 



\m/2\ 



[m- ]) 



j=0 



(j!)2(m-2j)! ■ 



(2.22) 



It is conjectured [39, Conjecture 2.1] that 



TriCyl(m) 



m 



^m ~r 



E 

d\m] l<d<m 



(j){m/d) td 



(2.23) 



where td is the coefficient of z'^ in the expansion of (1 -|- 2; -|- z^Y, i.e., the central 
trinomial coefficient (given as sequence A002426 in [30]), and 0(a;) is Euler's 
totient function. 

The values of all these dimensions for m < 16 are displayed in Table 2 of Ref. [6]. 



^ We have recently rederived this formula using a different method [38] . 



2.3 Beraha—Kahane— Weiss theorem 

A central role in our work is played by a theorem on analytic functions due to 
Beraha, Kahane and Weiss [16-19] and generalized slightly by one of us [20]. The 
situation is as follows: Let D he a domain (connected open set) in the complex plane, 
and let ai, . . . , aM, Ai, . . . , Am {M > 2) be analytic functions on D, none of which is 
identically zero. For each integer n > 0, define 

M 

Uz) = Y,M^)h{zr. (2.24) 

fc=l 

We are interested in the zero sets 

Z{fn) = {zeD: fn{z) = 0} (2.25) 

and in particular in their limit sets as n ^ oo: 

liminf 2(/„) = {z E D: every neighborhood U 3 z has a nonempty intersection 

with all but finitely many of the sets Z{fn)} (2.26) 

limsup2(/„) = {z E D: every neighborhood U 3 z has a nonempty intersection 

with infinitely many of the sets Z{fn)} (2.27) 

Let us call an index k dominant at z if |Afc(2;)| > |Ai(2;)| for all / (1 < / < M); and let 
us write 

-Dfc = {z E D: k is dominant at z} . (2.28) 

Then the limiting zero sets can be completely characterized as follows: 

Theorem 2.2 [16—20] Let D be a domain in C, and let ai, . . . , aM, Ai, . . . , Am 
(M > 2) he analytic functions on D, none of which is identically zero. Let us further 
assume a "no-degenerate-dominance" condition: there do not exist indices k ^ k' 
such that \k = ooXk' for some constant u with \lj\ = 1 and such that Dk (= D^') has 
nonempty interior. For each integer n > 0, define fn by 

M 

fc=l 
Then liminf 2^(/„) = limsup2(/„), and a point z lies in this set if and only if either 

(a) There is a unique dominant index k at z, and ak{z) = 0; or 

(b) There are two or more dominant indices at z. 

Note that case (a) consists of isolated points in D, while case (b) consists of curves 
(plus possibly isolated points where all the Afc vanish simultaneously). Henceforth we 
shall denote by B the locus of points satisfying condition (b). 

We shall often refer to the functions A^ as "eigenvalues" , because that is how they 
arise in the transfer-matrix formalism. 

10 



3 Structural properties of the square- lattice trans- 
fer matrix 

In this section we prove some structural results concerning the transfer matrices 
of square-lattice Potts models (and in particular chromatic polynomials) with free 
or cylindrical boundary conditions. We begin by proving some general results (Sec- 
tion 3.1) concerning the polynomial dependence in q of the transfer- matrix entries and 
the large-g behavior of the eigenvalues. Then we provide explicit closed-form expres- 
sions for the dominant diagonal entry of the transfer matrix with free or cylindrical 
boundary conditions (Sections 3.2 and 3.3). 

3.1 General results 

Let us begin by considering the case of the full Potts-model partition function. 
Indeed, we can be quite a bit more general, and consider any transfer matrix built 
out of operators of the form 

H = ^ CA]^Ji,i+i (3.1a) 

AC{l,...,m} i&A 

V = Y. ^^W^i (3.1b) 

BC{l,...,m} ieB 

with arbitrary coefficients {ca} and {ds}- (Recall that site m + 1 is to be understood 
as a synonym for site 1.) This general form includes as particular cases the square- 
lattice transfer matrix with free or cylindrical boundary conditions and arbitrary 
couphngs {vi^i+i} and {vi}. 

Let us now recall that "id" denotes the partition of {1, ... , m} in which each ele- 
ment is a singleton (i.e., {{1}, {2}, . . . , {m}}). Let us call a partition P of {1, ... , m} 
non-trivial if it is not "id". 

Proposition 3.1 For any operators H andV of the form (3.1), where the coefficients 
{ca\ and {ds} are numbers (i.e., independent of q), the diagonal entry ti^{m) of the 
transfer matrix T{m) = VH associated to the basis element Bid is a polynomial in q 
of degree at most m, of the form 

tid = C0(i{i^,,,^m} g™" -|- terms of order at most q"^^^ . (3.2) 

All other entries of the transfer matrix T{m) are polynomials in q of degree at most 
m — 1. 

Remark. If H is a projection, then the diagonal entry t[^{m) of the modified 
transfer matrix T'{m) = HVH associated to the basis element fid = Heid is given 
by tJd = tid, indeed, all the entries of T'(m) are equal to the corresponding entries 
of T(m). This follows immediately from the fact that Te-p = ^p, t-pp/ep/ implies 
T'fp = ^^, t-p-pii-pi. 

11 



Proof of Proposition 3.1. First of all, it is obvious that each entry in the 
transfer matrix T(m) is a polynomial in q. Indeed, from (2.7)/(2.9) it is clear that 
we get a factor of q every time we apply the operator Dj to a partition in which i 
is a singleton. We can maximize the number of factors of q by applying the vertical 
transfer matrix V to the vector eid that corresponds to the partition in which every 
site is a singleton. In particular, from (3.1b) we have 

Veid = ^ dsq^^^ Gid (3.3a) 

BC{l,...,m} 

= ('^{i,...,m}Q'™' + terms of order at most q^~^^ eid • (3.3b) 

If we apply the vertical transfer matrix to any other partition, we get a polynomial 
in q of degree at most m — 1. 

Let us now consider the quantity Heid: 



Heid = y^ ca TT Ji,i+i eid (3.4a) 




= C0eid+ 2_^ a-pe-p , (3.4b) 

■p non-trivial 

for some quantities a-p that are polynomials in {ca} (and of course independent of q). 
Using (3.3)/(3.4) it is obvious that 

VHeid = C(s ^ rfijgl^leid + ^a^(g)ep , (3.5) 

BC{l,...,m} V 

where the coefficients a'-p{q) are polynomials in q of degree at most m — 1. D 

In view of Proposition 3.1, we shall henceforth refer to tid (or t[^) as the "dom- 
inant diagonal entry" in the transfer matrix, as it is indeed dominant at large |g|. 
Furthermore, we can deduce from Proposition 3.1 the leading large-|g| behavior of 
the eigenvalues. We begin with a simple perturbation lemma: 

Lemma 3.2 Consider an N x N matrix M{C,) = {Mij{^))f-^ whose entries are 
analytic functions of C, in some disc \C,\ < R- Suppose that Mn = 1 and that Mij = 
O{0 foi" {hi) 7^ (1)1)- Then, in some disc |^| < R' , M(^) has a simple leading 
eigenvalue /i*(0 that is given by a convergent expansion 

oo 

/i,(0 = l + ^a^e' (3.6) 

k=2 

[note that cti = 0/ with associated eigenvector 

oo 

MO = e, + J2^,e, (3.7) 

fc=l 

while all other eigenvalues are 0{^). 

12 



The key fact here is that the eigenvalue shift in (3.6) begins at order ^^, not order 
Proof. We have 



N 



det[/x - M(0] = (/i - 1) ll[^^ - MuiO] + e'i^(/i,0 , (3.8) 



j=2 



where F{fi,^) is a polynomial in /i whose coefficients are analytic functions of ^ for 
|.^| < R. The polynomial P{fi) = (/i— 1) ni=2[/^~^««(0] has, for all sufficiently small 
|.^|, a simple root at ;U = 1 and roots (not necessarily simple) at /i = Mjj(^) = 0{^). 
The simple root at yU = 1 moves analytically [41] under the perturbation ^^F(/x, ^) — 
let us call this root fi*{C,) — and so is given by the convergent expansion (3.6). The 
corresponding eigenvector also moves analytically under the perturbation M{C,) = 
diag(l, 0, . . . , 0) + O(^), which proves (3.7). To see that all other eigenvalues are 
0{^), it suffices to consider the reduced matrix 

and observe that all its entries are 0{^). D 

Remark. The "small" eigenvalues need not be analytic in ^. For instance, 



M(0 = H (3.10) 




has eigenvalues /i = 1 and fx = ±,^^/^. D 

To apply Lemma 3.2 to our transfer matrix T, we set ^ = q~^ and M = T/tid. We 
then have: 

Corollary 3.3 Consider operators H and V of the form (3.1) where the coefficients 
{ca} and {ds} are numbers (i.e., independent of q) and C0(i{i^...^m} ^ 0. Then T has a 
simple eigenvalue K that is analytic for large \q\ and behaves there like C0(i{i^...^j„}g™.- 
more precisely, it has a convergent expansion 

, oo 

-^ = l + J2a,q~' (3.11) 



^id fc=2 



m— l^i 



[so that, in particular, A^ — tjd = 0(g™ ^)y. All other eigenvalues are 0{q 

Let us now return to the case of main interest, in which H and V are the transfer 
matrices (2.8)/(2.9) for the chromatic polynomial Vi = Vij = —1. In this case we can 
sharpen (3.2) by providing explicit expressions for the lower-order terms. We must 
now distinguish between free and cylindrical boundary conditions, and we shall treat 
each case in a separate subsection. 

13 



3.2 Free boundary conditions 

Let us consider a square-lattice grid of fixed width m > 1 and free boundary 
conditions. Let us also assume that all horizontal edges have weight v and all vertical 
edges have weight v'; they need not be —1. The horizontal transfer matrix (2.8a) is 
thus 

m—l 

H = J](l + t;J,,,+i) (3.12) 

i=l 

and the vertical transfer matrix is 

m 

V = l[{v'I+Di). (3.13) 

1=1 

Consider first the action of H on the start vector eid. It generates 2™~^ terms, 
each of which corresponds to a partition V in which all the blocks are sequential sets 
of vertices in {l,...,m} (we shall call these sets "polymers"). Furthermore, each 
polymer of size i picks up a factor v^~^. 

Consider next the action of V on a basis vector e-p corresponding to an arbitrary 
partition V = {Pi, . . . , Pk}- If we are to end up with the partition Bid, then for each 
block Pj we must either choose the delete operator Dj for all i G Pj (the last deletion 
gives a factor q) or else choose the delete operator for all but one i G Pj and choose 
v'l for the last site (this can be done in \Pj\ ways). We therefore have 



Ve 



r 



k 

b=i 



Bid + other terms (3-14) 



where "other terms" means terms involving ep/ with V' ^ id. We thus obtain a factor 
q + iv' for each block of size i = \Pj\ > 1. 

Putting these facts together, we conclude that 

VHeid = tF(m)eid + other terms , (3.15) 

where tFlm) is the partition function for a one-dimensional m-site polymer gas (with 
free boundary conditions) in which each site must be occupied by exactly one polymer, 
and each polymer of length i > 1 gets a fugacity fii = v^~^{q + iv'), i.e. 

oo k 

tF{m) = J2 E U'^'^-^q + hv') . (3.16) 

il + . . . + £k = m 

To solve this polymer model, let us introduce the generating function ("grand parti- 
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tion function" 



<^v(z) 



Y,z-t 



F[.m) 



m=l 



fc=i £i,...,4>ii=i 



1 - ^f^) 



where 



^{z) = ^z'v'~\q + lv') = z 



e=i 



1 — zv 



:i - zv) 



is the total weight for a single polymer of arbitrary size. We therefore have 

^ , , (q + v')z - qvz'^ 

^■p(z) = — 

^' l-{q + 2v + v')z + v{q + v)z^ ' 

When f = f' = — 1 this reduces to 



(3.17a) 

(3.17b) 
(3.17c) 

(3.18) 

(3.19) 
(3.20) 



l-{q-3)z-{q-l)z^ ' 

As a check we have expanded (3.20) in powers of z, and verified that it agrees with 
available dominant diagonal elements tF{m) for 1 < -m < 12 (see [6] for m < 8). 

The next step is to get an explicit expression for tp{m). Using the notation 
[z"^]P{z) for the m-th coefficient in a polynomial or power series, we have the alternate 
representations 



1 — az — hz"^ 



Lm/2J 

E 

j=0 



Y, r'^\a^-'^v 



(3.21a) 



2-™ Y. C". ) «'""''■(«' + 46)-'' . (3.21b) 



The first of these comes directly from [2;™] X]fclo('^'^ '^ bz'^)'', while the second is ob- 
tained by factoring the quadratic and using partial fractions. Using (3.20) we have 

*,(,„) = |."H.,(.) = (,-l)|--.|_l_5 + ,|,...-|_1_ (3.22) 
where a = q — 3 and h = q — 1. Inserting this into (3.21a) we have 

L(m-1)/2J . _ 1 _ -X 

j=o V J / 



L(m-2)/2j 

E/m — 2 — 7 
( J 

15 



(?-3) 



m-2-2jr' 



(g-1)^ 



(3.23) 



which is manifestly a polynomial in q of degree at most m. Furthermore, the term 
of order g™ comes only from j = in the first sum and has coefficient 1 (thereby 
confirming explicitly what we already knew from Proposition 3.1), so the degree is 
exactly m. 

We may therefore write tpiTn) explicitly as a polynomial in q with certain coeffi- 
cients a^{m): 

m 

tpM = J2{-l)'^al{m)q^-K (3.24) 

fc=0 

The next step is to obtain a closed formula for the coefficients alim) with m > 1 
and < k < TTi. We shall prove that, for each fixed A; > 0, the coefficient alim) is in 
fact a polynomial in m of degree exactly kJ We begin by expanding the binomials in 
(3.23): 

al{m) = (-l)'^[g™-'^]tF(m) (3.25a) 

-(-!)' "I" g (" ■; " ') (::i."-';) C ; ') (-3)— (-d--' 
- (-)' "T t C" -; - ') ( j_?_-/i (-3)— (-1)- 

(3.25b) 

L(m-2)/2j «,/ oo-\/- 

m — 2 — j\ f m — 2 — 2j \ f j 

(3.25c) 



i=o i=o 



E 5:('"7"-'j(;:,_r-\)(;)(-i''^'"-'"' 



Let us consider the first sum in (3.25c): 

[(m-l)/2j oo 

= E E4'^^(^'^)- (3.26) 

The goal is to substitute X^jl^^ by something independent of m, e.g., X]j=o- 
Indeed, if A; = [(m — 1)/2J, the identity is trivial. Let us consider next the case 
k < \_{m — l)/2j. Then, if we prove that the sum 

L(m-1)/2J 

J2 Sf}{m,k) = 0, (3.27) 

j=k+i 



More precisely, a\{ra) is the restriction to integers m > max(fc, 1) of such a polynomial. 
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then, we have not modified the resuh of (3.26) by changing the upper index in the 
sum over the variable j. The second binomial appearing in (3.26) vanishes whenever 
k + i— 1 — 2j<0. On the other hand, the third binomial is non-vanishing only if 
j > i-l. Therefore, if j > k and j > i-1, we have k+i-l-2j < k+i-l-k-{i-l) = 
0. So, S^^\m,k) = whenever j > k.^ Finally, let us consider the third case 
k > \_{m — l)/2j. Then, by making this change in the upper index in the sum over j, 
we are adding some extra terms 

k 

E 4VK^)- (3-28) 

i=L(m-l)/2j+l 

In this case we have to focus on the first binomial of (3.26). This binomial is nonzero 
only when < j < \_{m — l)/2j or when j > m. The first of these do not appear in 
(3.28); and since m > k, the second appears only when j = k = m. The contribution 
of this extra term is 1. Thus, (3.26) reduces to 



j=0 £=0 



772 — 1 — A / 772 — 1 — 2j \ f j + 1 

J )\k + ^-l-2j 



km 



Let us now consider the second sum in (3.25c): 



(3.29) 



[(m-2)/2j oo 



m — 2 — j\ f m — 2 — 2j \ fj 



[(m-2)/2j oo 

= E E^]?(^'^)- (3-30) 

The goal is now to substitute X]i'=o^ ^^ Sj=o- ^^ before, the case k — 1 = 
[(m — 2)/2j is trivial. Let us now suppose that k — 1 < \{m — 2)/2j. Then, the 
second binomial vanishes whenever k + d, — 1 — 2j < 0, the third binomial is non- 
vanishing only if j > i. Therefore for j > k — 1 and j > i we have k -\- i — 1 — 2j < 
k-\-£ — 1 — {k — 1) — i = 0. So S-l{m, k) = whenever j > k — 1. Let us finally 
consider the third case k — 1 > [{m — 2)/2j. Then, we should consider the extra 

terms 

fc-i 

E S^^{m,k). (3.31) 

i=L(m-2)/2j+l 

In this case we have to focus on the first binomial of (3.30). This binomial is nonzero 
only when < j < [{m — 2)/2j or when j > m — 1. The first of these do not appear 
in (3.31); and since m > k, the second appears only when j = k — l=m — 1. The 



^This is true even if m is treated as an algebraic indeterminate. 
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contribution of this extra term is again 1. Thus, (3.30) reduces to 

(3.32) 
Putting together (3.29)/(3.32), we find that the two contributions 6km cancel 
exactly, and that (3.25c) can be written as 

(3.33) 

where the independent variable m does not appear in the summation limits. After 
some straightforward but lengthy algebra we can rewrite the above formula in the 
more compact form 

From (3.33) or (3.34) we see that a^im) is (the restriction of) a polynomial in m 
of degree at most k, as m appears only in the upper index of the binomials and 

"'rn^ ml m{m — l){m — 2) ■ ■ ■ {m — j + 1) 

( O.OO ) 



J/ J! Ji 

is a polynomial in m of degree j. [Here and in what follows, we use Knuth's [42] 
notation for falling powers: xi = x{x — l){x — 2) ■ ■ ■ (x — j + 1) .] 

The degree of the polynomial al{Tn) is in fact exactly k. To see this, let us extract 
the term of order m'' from (3.34). The second sum in (3.34) does not contribute, as 
it is a polynomial in m of order at most k — 1; the only contribution comes from the 
first sum: 

^ / 1 \pok—p nk 

We can summarize all this into the following proposition: 

Proposition 3.4 Let H and V be the transfer matrices (2.8a)/(2.9) for the chro- 
matic polynomial Vi = fi,i+i = —1 with free boundary conditions. Then the dominant 
diagonal entry in the transfer matrix can be written as 

m 

t^m) = $^(-l)^a^(m)g™-^ (3.37) 

fc=0 
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where each a^{Tn) is a polynomial in m of degree k given by (3.34)- 
For the first few values of k, we obtain 



aQ[m 



F/ 

a^ [m 

F/ 
F/ 
F/ 
F/ 

05 [m 

F/ 



1 (3.38a) 

2m - 1 (3.38b) 

2m^ - 3m + 1 (3.38c) 

|m3 - 4m2 + |m (3.38d) 

2^4 _ w^3 ^ ^^2 ^ 5^ _ 1 (3.38e) 

^m^ - 2m* + ^m^ + fm^ - ^m - 2 (3.38f) 
^m^ - i|m5 + f m^ + Jm^ - ^m^ - ||m + 12 (3.38g) 



We can also prove the following result concerning the polynomials al{m):^ 
Proposition 3.5 For each integer s > 1, the quantity 



al{k - s) 



p=0 \ f / j.=0 ^ 



p+l 

k — p — r 



p=i ^ 



k — s — 1 — p\ Y^^j. 



Y,3' 



k — s — 2p\ f p — 1 



p — 1 / -^^ V r I \k — p — r 



(3.39) 



is, when restricted to k > s + 1, given by a polynomial in k of degree max(0, s — 3), 
with leading coefficient 



[k'-']alik-s) 
and first subleading coefficient 

r->i{k - s) 



\s+l 



(5-3)! 



-lYs 



for s > 3 



(3.40) 



for s > 4 . 



(3.41) 



2(s-4)! 
The proof of this Proposition can be found in Appendix A. 2. For the first values 



^This Proposition refers, however, to the regime m < k that does not contribute to the sum 
(3.37). 
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of s, we have 






al{k-l) = 


for A; > 2 


(3.42a) 


alik-2) = 


for A; > 3 


(3.42b) 


alik-3) = 


1 for A; > 4 


(3.42c) 


al{k-4) = 


-{k - 2) for A; > 5 


(3.42d) 


al{k-5) = 


i(A;2-5A;-2) for A; > 6 


(3.42e) 


al{k-6) = 


-|(A;3-9A;2-4A; + 66) 




= 


-|(A;-3)(A;2-6A;-22) for A; > 7 


(3.42f) 


al{k-7) = 


2^(A;4-14A;3-A;2 + 350A;-384) for A; > 8 


(3.42g) 


al{k-8) = 


^l{k^ 20A;^ + 15A;3 + 1100A;2 3016A; 2400) 





= Yio('^-4)(A;4-16A;3-49A;2 + 904A; + 600) for A; > 9 (3.42h) 

In particular, we see from (3.42a,b) that (m — A; + l)(m — A; + 2) is a factor of af (m) 
for A; > 3. 

We find empirically that, for each integer s > 3, the polynomial p^{k) that matches 
(— l)*^"'^(s — 3)!af(A; — s) for integer k > s + 1 has all integer coefficients; and we 
further find empirically that for even integers s > 4, we have p^{s/2) = 0, so that 
the polynomial p^{k) has k — s/2 as a factor. 

We also find 



al{k) 
al{k + l) 



= 


F2k for A; > 1 


(3.43a) 




(2A; + l)F2,+2 - (A: - 4)F2,+i ,^ ^^^, 




5 




1 

/5 


V 2 J [ 2 ) 


(3.44) 



where 



are the Fibonacci numbers: see [30, sequences A001906/A088305 and A038731]. We 
prove (3.43a) in Appendix A. 4. We have checked (3.43b) up to A; = 100, but do not 
have any proof. ^° Please note that (3.43a,b) give the low-order coefficients in the 
polynomials tF^m): 



al{k + l) 



:-mq'Mk+i] 



(3.45a) 
(3.45b) 



Since we have proven that of (m) is a polynomial in m of degree k, it is also of 
interest to obtain explicit expressions for the coefficients of this polynomial, which we 



^°With some more work it might be possible to find a proof of (3.43b) using the same strategy as 
was used for the proof of (3.43a) in Appendix A. 4. 
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write as 



al{m) 



^ ( 1^^oA:-2£+l 



£=0 



(A; -£)!(£ + 2)! 



(3.46) 



here the prefactors have been chosen to make many (though not all) of the coefficients 
a\^ integers (in fact, all of them are integers for £ < 5, see below). Now we use the 
well-known expansion of the falling powers in terms of Stirling cycle numbers [42], 



X- 



E 

c>0 



{-ly-^x^ , 



(3.47) 



and expand all the binomials in (3.34) involving m. We arrive after some algebra at 
the following expression: 



(^kl — 



/_lY2k-2e+i 
(k -£)]{£ + 2)\{-l)'' 



m''-%l(m) 



(3.48a) 



2fe-2£+l 



Z-^ Z-^ n r- ^-^ 



a 

a=0 c=0 >- -^ 
k k—p 



k k—p 

EE 

p=0 r=0 

k-t 



p+l \ (-3)' 
k — p — rj p\r\ 



'- -> d=0 



k — i — dj \d 



'i+2pr-\i+p) 



a+d—k-\ 



/-^ Z-^ \ h — n — r (n — ^\\r\ ^-^ ^-^ 



p — 1 



[sy 



^k — p — rj (p — l)!r! 

p=l r=0 ^ ^ ^ ^-^ ' a=0 c=0 

k-e 

E 

d=0 



p — 1 

a 



X 



k 



-d 



{2p) 



c—d/ 



p) 



a+d-k+e 



(3.48b) 



By computing (3.48b) for integers A; > £ > 0, we find empirically that a^^ is in fact, 
for each fixed i, (the restriction of) a polynomial in k of degree i. The ffist few of 
these polynomials are: 



F 
O'kfl 


1—1 






(3.49a) 


F 
«fc,l 


= 3A; + 3 






(3.49b) 


F 
«fe,2 


= 6k^-Uk + 52 






(3.49c) 


F 

afe,3 


= lOA;^ - lOOA;^ + 130A; + 240 






(3.49d) 


F 
«fe,4 


= I5k^ - 330k^ + 845P - 18A; - 1928 






(3.49e) 


F 
«fc,5 


= 2lk^ - 805k^ + 5005P + 749A;2 + 8358A; - 


87360 


(3.49f) 


F 

«fc,6 


- 28A;^ 1652P + 20020A;^ ^^^^^^e- 


1 278096 j^2 
o 






1 3141872 u 3838336 
1 9 ''' 3 






(3.49g) 


The fact that a^g 


= 1 for all A; > is just a restatement of 


(3.36) [compare (3.46) . 
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3.3 Cylindrical boundary conditions 

Let us now consider a square-lattice grid of fixed width m > 1 and cylindrical 
boundary conditions. (Please note that for m = 1 the horizontal edges are loops, and 
that for m = 2 there are two horizontal edges connecting the pair of sites in each 
row.) Let us also assume that all horizontal edges have weights v and all vertical 
edges have weights v'; they need not be — L We proceed analogously to the preced- 
ing subsection, making the changes necessary to handle cylindrical rather than free 
boundary conditions. 

Consider first the action of H on the start vector Bid. It generates 2™" terms, 
each of which corresponds to a partition V in which all the blocks are sequential sets 
of vertices on the m-cycle (we shall call these sets "polymers"). Furthermore, each 
polymer of size i < m picks up a factor f^^^, while a polymer of size m picks up a 
factor v"^ + mv'^~'^ (the v^ comes from the case in which all edges are occupied, while 
the mv"^~^ comes from the m cases in which all edges but one are occupied). The 
action of V is identical to that for free boundary conditions. 

The upshot is that we have 

VHeid = tp(m)eid + other terms , (3.50) 

where tp(m) is the partition function for a polymer gas on the ?T2-cycle in which each 
polymer of length i > 1 gets a fugacity 

^ _ jv^-\q + iv') ioTl<i<m-l 

I v'^~'^{v + m)(q + mv') for £ = m 



Please note that 



/2, = <(^^ forl<£<m-l ^3^^2) 

(t> + m)firn ioT i = m 



where /i^ are the fugacities for free boundary conditions considered in the preceding 
subsection. 

We can obtain the tp(m) by using a simple recursion relating the periodic and 
free cases: 

m—l 

tp{m) = y^kfik'tpi'm — k) + Jim ■ (3.53) 

fc=i 

To see this, single out a site (e.g. 1) and let k > 1 be the size of the polymer placed 
on it. li k < m — l, we have k ways of placing this polymer such that the selected site 
belongs to it, with fugacity /i^ for each such placement; and for the rest of the ring, 
the total weight of all admissible polymer configurations is simply t-F{m — k). Finally, 
if /c = 772, there is only one way of placing the polymer, and it receives fugacity /im. 
This proves (3.53). 

In order to compute explicitly the tp(m), it is convenient to introduce the gener- 
ating function 

oo 

$p(2) = ^z™tp(m) . (3.54) 

m=l 
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Note next that the upper hmit on the sum in (3.53) can be changed to cx), provided 
that we define t^ii) = for i < [which is anyway imphcit in the definition (3.17a) of 
the generating function <I>f(-2)]- Multiplying both sides of (3.53) by z™ and summing 
over m, we arrive easily at the equation 

$p(z) = z—^[1 + <^f{z)]+v'^{z) (3.55a) 

z d'^iz) , ^ / , X 

+ v'^{z) (3.55b) 



1 - ^(2) dz 

where "^{z) = Xl^i -^^/^^ i^ defined in (3.18)."'^^ When v = v' = —1, we obtain the 
final formula 

^^ \l + zj l-{q-3)z-{q-l)z^ ^ ' 

_ qz^ + {q + l)z + 2 2-{q-3)z 

{1 + zy l-{q- 3)z ~{q- \)z^ ' ^^'^^^^ 

By expanding this function in powers of 2;, we have checked that it agrees with the 
known dominant diagonal elements tp{m) for m < 13 [6,7]. 

It is now easy to extract the partition function tp(m): using (3.21a) we get 

tp{m) = [z"']<^p{z) = (-l)™(g-m-2) 

[m/2j 



j=0 \ J / 

- E r ' (^-3r-^^■(g-l)^ (3.57) 



which is manifestly a polynomial in q of degree m. 

We can now define the coefficients a^{iTi) in the same way as for free boundary 
conditions: 

m 

tp(m) = J2{-lfal{m)q"^-\ (3.58) 

fe=0 

where k and m are integers satisfying m > 1 and < k < m. However, it is slightly 
more convenient to extract explicitly part of the term (— l)'"(g — m — 2) from (3.57), 
and define d^im) to be the coefficients in what remains: 

m 

tp{m) = {-ir{q - m - 1) + ^(-l)^^ Q"'''' ■ (3-59) 

A;=0 



^^Note that the term of order z^ in (3.55) vanishes whenever v = —1 (irrespective of the values of 
q and v'). This reflects the fact that fp(l) = whenever v = —1 because of the loops at each vertex. 
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Notice that the relation between a^im) and a^im) is rather simple: for fixed m > 1 
we have that 

{d^ (m) for < A; < m — 2 

d^{m) — 1 for A; = m — 1 . (3.60) 

d^{m) — (m + 1) for k = m 

Expanding the binomials in (3.57), we have 

-TEr-;-0G:7-^g(;)M)'3'-. 

- 5km (3.61) 

Again we want to substitute the m-dependent upper index in the sum over j by 
something independent of m: e.g., by k. 

In the first sum there are two non-trivial cases: a) If A; < \m/2\ , then the second 
binomial vanishes whenever k + i — 2j < 0, and the third binomial is non- vanishing 
only if j > i. Therefore for j > k and j > i we have that k + i — 2j < k + i — k — i = 0. 
So all these terms vanish, b) If A; > \m/2\ , then the first binomial does not vanish 
only when < j < \m/2\ or when j > m + 1. As we are adding terms with 
[m/2j +1 < j < k — 1 < m—1, none of them give rise to a non- vanishing contribution. 

In the second sum we play a similar game: a) If A; < [(m — 1)/2J, the binomials 
involved are the same as for the first sum, so the same result applies here too. b) 
If A; > \_{m — l)/2j, then we are adding terms with \_{m — 1)/2J + 1 < j < A; — 1. 
The first binomial does not vanish only when < j < [{m — l)/2j or when j > m. 
The first of these do not appear in the extra terms; and since m > k, the second 
appears only when j = k = m, giving rise to an extra contribution equal to 1. This 
contribution cancels out exactly the term —6km in (3.61). 

Putting all the pieces together, we end with the following expression for d^{m): 



k 



3=0 



m — j\ fm — I — j 
J J \ J 



a,v)^E(-i)' ^r7]-r~r') e ( ,: ^_'y ;>"'- 



(3.62) 
where the independent variable m does not appear in the summation limits. After 
some straightforward but lengthy algebra we can rewrite the above formula in a more 
compact form: 



Kr)^i:M)'?r;!rOs 



k — p — r^ 

(3.63) 
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It is clear from (3.63) that a^{fn) is (the restriction of) a polynomial in m of degree 
at most k. To see that its degree is exactly k, let us extract the term of order m'': 



m'Kim) 



3^ A {-1)P3^-P 



k\ 



3fc 2^ - 3^ 
I\^ k\ 



k\ 



+ 0. 



(3.64) 



Let us also remark that the constant term in a^^,{v(l) vanishes whenever /c > 1: 



m ]%(m) = 4o 







if A; = 
if A;> 1 



(3.65) 



We can summarize the foregoing results in the following proposition: 



Proposition 3.6 Let H and\l he the transfer matrices (2.8h)/(2.9) for the chromatic 
polynomial Vi = fj,i+i = — 1 with cylindrical boundary conditions. Then the dominant 
diagonal entry in the transfer matrix can be written as 



tAm) = i-iriq-m-l)+J2i-^Mirn)q"^-'' 



(3.66) 



fc=0 



where each a^{m) is a polynomial in m of degree k given by (3.63). 
The first polynomials a\{m) are given by 



Cq [m 

~P/ 
a^ [m 

-P/ 

~P/ 
Og [m 

~v / 
a^ [m 

-P/ 
Og [m 

~P/ 
Og [m 



1 

2m 

2rn^ — m 



|m^ — 2m^ — |m 

|m^ — 2m^ — |m^ + |?t2 

Y^m^ — |m*^ + |m^ + ym^ — |m 

45 3 9 2 45 6 



We can also prove the following result concerning the polynomials ai,{m 
Proposition 3.7 For each integer s > 1, the quantity 



(3.67a) 
(3.67b) 
(3.67c) 
(3.67d) 
(3.67e) 
(3.67f) 
(3.67g) 

1.12 



a^ik-s) 



k — s\ s-^, _,,„ k — s fk ~ s — p — l\ sr-^ ^fk — s — 2p 



k 



^ p=i 



p 



p — 1 



k~p 
r=0 



p 

k — p — r 
(3.68) 



^^This Proposition refers, however, to the regime m < k that does not contribute to the sum 
(3.66). 
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is, when restricted to k > s, given by a polynomial in k of degree s, with leading 
coefficient 

[kral{k-s) = tip (3.69) 

and first suhleading coefficient 

r~']al{k-s) = ^^^^^- (3.70) 

Furthermore, 5.^(0) = 0, so that the polynomial representing a^{k — s) for k > s has 
a factor k — s. 

The proof of this Proposition can be found in Appendix A. 3. For the first values 
of s, we have 

al{k-l) = k-1 ioTk>l (3.71a) 

dl{k-2) = -l(A;2-3A; + 2) 

= -|(A;-l)(A;-2) for A; > 2 (3.71b) 

al{k-3) = li^ - 6k^ + 5k + 12) 

= l{k + l){k - 3){k - 4) for A; > 3 (3.71c) 

al{k-4) = -^{k^-We + lie + Mk-lQS) 

= -j^{k + 3){k-2){k-4:){k-7) for A; > 4 (3.71d) 

hl{k-5) = Y^(A;5-15A;4 + 25A;3 + 375A;2-1346A; + 480) 

= ^{k-5){k'^-10k^ -25k'^ + 250k-96) for A; > 5 (3.71e) 

al{k-6) = -yiQ(A;6-21A;5 + 55A;^ + 1065A;3-6176A;2 + 3636A;+ 15840) 

= -^(A;-3)(A;-6)(A;4-12A;3-71A;2 + 642A; + 880) for A; > 6 

(3.71f) 

We find empirically that, for each integer s > 1, the polynomial p^{k) that matches 
(— l)*~^s! a^{k—s) for k > s has all integer coefficients; and we further find empirically 
that for even integers s > 2, we have pf (s/2) = 0, so that the polynomial p^{k) has 
k — s/2 as a factor. 

We also find empirically 

al{k) = F2k+i + F2k-i - I (3.72a) 

alik + l) = E T 1^^^ (3-72b) 

i=o V -^ / 

where F„ are again the Fibonacci numbers (3.44): see [30, sequences A005592 and 
Al 17202]. We have checked these relationships up to A; = 100, but do not have any 
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proof. ^^ Please note that (3.72a,b) give the low-order coefficients in the polynomials 
tp(m): 



al{k) = (-l)1g%(A;) + (A; + l) 
~al{k + l) = (-l)'=[gi]tp(A; + l) + 1 



(3.73a) 
(3.73b) 



Since a^{m) is a polynomial in m of degree k, we are again interested in obtaining 
the coefficients a^^ defined by 



al(m) 






-1Y2 



er)k-2e+i 



{k-i)\{i + 2)\ 



^lem^ ^ 



(3.74) 



With the help of (3.47), we obtain after some algebra the following result: 



~p 



2^-2^+1 {li 

k k—p j 



k 
k~ 

-3)" 



p— 1 r 

TrT^ \ k — T) — T J ^ ^ 

p=l r=0 ^' • \ ^ / a=0 c=0 



k-e-i 

d=0 



P 



d \k-i-d- 1 



p — 1 
a 



{2py~''{i+p) 



a-k+e+d+1 



(3.75) 



By computing (3.75) for integers k > i > 0, we find empirically that a^^ is in fact, 
for each fixed i, (the restriction of) a polynomial in k of degree i. The first few are: 



~p 
~p 

~P 

'^k,2 

~P 
^k,3 

~P 
'^k,A 

~P 
~P 



1 

3A;-3 

6A;2 - 38A; + 52 

lOA;^ - 160A;2 + 390A; + 

15A;^ - 450A;3 + 2405A;2 - 1450A; - 7688 



(3.76a) 
(3.76b) 
(3.76c) 
(3.76d) 
(3.76e) 



2lk^ - 1015A;^ + 10465P - 1612lP - 37030A; - 151200 (3.76f) 



28A;^ - 1988P + 34580A;^ - mfl^k^ _ M^Zep 



480752 
9 



k~ 



1902976 



(3.76g) 



The fact that a^Q = 1 for all A; > is just a restatement of (3.64). 



^■^With some more work it might be possible to find a proof of (3.72a,b) using the same strategy 
as was used for the proof of (3.43a) in Appendix A. 4. 
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4 Large-g expansion of the leading eigenvalue 

In this section we compute the large-g expansion of the leading eigenvalue Av,(m) 
for both free and cylindrical boundary conditions, and determine empirically some of 
its remarkable properties. In Section 5 we shall provide theoretical explanations of 
some (but not all!) of these empirical observations. 

4.1 Overview of method and results 

In the preceding section we computed in closed form the dominant diagonal entry 
in the transfer matrix, tid, for a strip of width m > 1 with either free or cylindrical 
boundary conditions (denoted tp and tp, respectively). We found that this entry is 
in each case a polynomial in q of degree m: 

m 

t^m) = J2{-l)'al{m) q--' (4.1a) 

k=0 



tp(m) = J2^-l)'al{m)q^-' (4.1b) 

k=0 

m 

= (_i)™(^ - m - 1) + J2i-^f~^lM r-' (4.1c) 

fe=0 

We furthermore computed in closed form the coefficients a^(m) and a^{m), which are 
in fact polynomials in m of degree k [cf. (3.34) and (3.63)]. For instance, the leading 
few terms for large |g| are 

tp{m) = g" - (2m - l)g™-^ + {2m^ - 3m + l)q'^-^ 

- (|m^ - Am^ + ^m)q'^-^ + ... [for m > 1] (4.2a) 
tp{m) = g" - 2mg'"-i + (2m2 - m)g™-2 

— (|m^ — 2m^ — I 



-m^ — 2m^ — ^m)q''^' ^ + . . . [for m > 5] (4.2b) 



In this section we want to carry out an analogous computation for the dominant 

P/P 

eigenvalue of the transfer matrix, which we call A* . Already from Corollary 3.3 we 
can conclude that A^(m) has, for large |g|, a convergent expansion in powers of g~\ 

A."(m) = 5^(-l)^&^(m)g-^ (4.3a) 

fc=0 

oo 

AfM = E(-l)'^^(^)^'""' (4.3b) 

fc=0 

and that the first two terms in this expansion coincide with those in the dominant 
diagonal entry: 

AF/P(^) _ tp/p(^) = 0{q"^-') (4.4) 
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and hence 

6^/P(m) = a'J^im) for A; = 0, 1 . (4.5) 

Here we shall go further and compute the coefficients bj (m) for I < m < 12f, 13p 
and < A; < 40.^^ Somewhat surprisingly, we shall find (for the m values we were 
able to study) that 



so that 



blim) = 


al{m) 


for k = 2,3 


(4.6a) 


blim) = 


al{m) 


for A; = 2, 3 and m > A; + 2 


(4.6b) 


Af(m) - 


- tpim) 


= 0(g-^) 


(4.7a) 


Ar(^) - 


- tp(m) 


= 0(g™-^) for m > 5 


(4.7b) 



rather than merely 0{q"^~'^) as Corollary 3.3 shows. We conjecture that this behavior 
holds for larger m as well. 

For some purposes it is slightly more convenient to use, in place of the coefficients 
b^{m), the modified coefficients b^{m) defined by 

oo 

Ar(m) = i-iriq - m - 1) + ^^(-l)'^^^^ 9"^"' (4-8) 

fc=0 

[analogously to (4.1c) for tp{m)]. Note that the relation between the coefficients 
b^{m) and b^{m) is the same as for the coefficients a^{m) and a^{m) [cf. (3.60)]. 

Most importantly, however, it is enlightening to pass from the eigenvalue Aj, (m) 
to its logarithm, which is a free energy, and define 



\og^^ = j:o'f{m)q-^. (4.9) 

^ k=i 



For cylindrical boundary conditions it is slightly more efficient to define the modified 
coefficients c^{m) by 

^ k=l 

In this section we shall see empirically that the coefficients Ck{m) behave in a 
much simpler way than the bk{m): namely, while bk{m) is, for large enough m, (the 
restriction of) a polynomial of degree k in m, we shall find that Ck{m) is, for large 



^^It would not be difficult to extend this computation to much larger values of A:, if we really 
cared. Extension to larger values of m is, however, an extremely demanding computational task. 
The dimension of the transfer matrix for m = 12f is 2947; for m = 13f it is 7889, which is beyond 
the capabilities of our current computer facilities. 
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enough m, (the restriction of) a polynomial of degree 1 in m. In Section 5 we shall 
discuss the theoretical interpretation of this empirical observation. 

We shall proceed as follows: Using the methods of [6, 7] we shall compute the 
transfer matrices for strips of width m < 12 for free boundary conditions and m < 13 
for cylindrical boundary conditions. ^^ From these we can extract the dominant eigen- 
value as a power series in q~^, i.e. for each available m we can easily compute as many 
coefficients bkim) and Ckim) as we please. "^^ We then observe empirically that, for 
each A; > 0, the coefficient bk{m) [resp. Ck{m)] is a polynomial Bj^ [resp. Ck] in m of 
degree k [resp. degree 1] provided that we restrict to integers m > some mmin(A;).^'^ 
Assuming that this empirical observation is accurate (i.e., that the polynomial be- 
havior persists to all larger m), we can infer the expressions for the polynomials Bk 
and Ck for A; < 31 (resp. k < 16) for free (resp. cylindrical) boundary conditions. 

4.2 Free boundary conditions 

Using the methods just described, we have obtained the leading eigenvalue A^(m) 
for < m < 12 as a power series in q~^ [cf. (4.3a)] through order k = 40. The resulting 
coefficients bl{m) are displayed in Table 1, and the corresponding coefficients c^{m) 
[cf. (4.9)] are displayed in Table 2. It is interesting to note that for all {k,m) that 
we have computed (i.e., 1 < -m < 12 and < A; < 40), the coefficients bl{m) and 
kc^{m) are integers. We observe empirically that, for each fixed k, the coefficients 
6^(m) are the restriction to integers m of a polynomial B^ in m of degree k, and that 
the coefficients c^{m) are the restriction to integers m of a polynomial Cf in m of 
degree 1, provided that we restrict attention to m > mf^:^j^{k) with 

if < A; < 6 
Kjk) = r.. • (4.11) 

- 2 if A; > 7 

Below this threshold ^«min(^)' ^^^ coefficients deviate from polynomial behavior. With 
our available data together with a few tricks described below, we are able to determine 
these polynomials for < A; < 33. 

First we start by trying to fit the coefficients b^{m) with m > iTL^^^^{k) to a 
polynomial Bj in m of degree k. As we need A; + 1 coefficients for such a polynomial. 




^'"'In fact, this was already done in ref. [6] for jti < 8 with both boundary conditions, and in ref. [7] 
for 9 < m < 13 with cyhndrical boundary conditions. Therefore, the only new transfer matrices we 
need to compute here are m = 9, 10, 11, 12 with free boundary conditions. See also Sections 5.2 and 
6 below for further results from this computation. 

^®To compute the dominant eigenvalue as a power series in q^^, we have applied the power 
method [43, Section 7.3.1] in symbolic form to the transfer matrix. Each iteration gives one additional 
term in the expansion of the dominant eigenvalue in powers of q~^. We can therefore compute the 
exact expansion up to any desired order in a finite number of steps. 

^^By contrast, for the dominant diagonal entry we have proven that a^(m) and d^{m) are poly- 
nomials in m of degree fc; and in this case the polynomial form holds for all allowable integers to, 
i.e. TO- > max(/c, 1). 
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we are able to obtain these polynomials only up to /c = 8. Please note that in all 
cases we have at least one data point more than the number of unknowns, so every 
fit can be tested at least on one extra data point. Our results are: 



Bl{m 
Bl{m 
Blim 
Bl{m 
B^(m 



B^{m 



Bl{m) 



2m- 1 

2m^ — 3m + 1 

,2 , 8 



4 S 



Am^ + |m 



2^4 _ W^3 + 23^2 ^ n^ _ 3 

3 3 D D 

^m^ — 2m^ + ym^ + |m^ — ^m + 9 



45 

8 
315 



m 



574 „ 



14^5 + 1^4 ^ ^^3 _ 2|3^2 ^ ^^ 



m 



15011 
210 



ifm^ + fm^ + f m^ 
m + 29 



1144 3 1 5947 ^2 
45 '"' "^ 90 '" 



315'" 



A^^ + 1^6 + 1^5 _ ^^4 + 4129^3 



302017 2 I 9041 _ aq. 
2520 '"' ^ 84 "' ^^ 



(4.12a) 
(4.12b) 
(4.12c) 
(4.12d) 
(4.12e) 
(4.12f) 
(4.12g) 

(4.12h) 

(4.12i) 



Notice that the three highest-order coefficients agree with those of the corresponding 
polynomial a\{m), i.e. 



Bl{m) 



af (m) for < A; < 3 

al{m) + 0{m^^'^) for A; > 4 



(4.13) 



However, there is a better way of extracting the desired information from our 
numerical data: instead of using the coefficients h\{m) as our basic quantities, we can 
use the related coefficients c\{rn) [cf. (4.9)]. The latter coefficients are empirically 
found to be, for each fixed fc, the restriction to integer ?72 of a polynomial Cf in m of 
degree 1, i.e. 

Cl{m) = alm + (3l, (4.14) 

provided that m > the same rn^^^{k) defined in (4.11). As we need only two coeffi- 
cients for such a polynomial (i.e., af and f3^), we are able to obtain these polynomials 
up to A; = 24 (if we want at least one extra data point to test the fit) or A; = 26 (if 
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we don't). Our results for k < 8 are: 

Cf{m 



C!i 



m 



Cf(m 



-2m + 1 



-m 



1 2 



|m- ^ 



f^-f 



/V-y^ ^O 

6 3 



^m ^ 



§lm — 243 
4 '"' 8 



(4.15a) 
(4.15b) 
(4.15c) 
(4.15d) 
(4.15e) 
(4.15f) 
(4.15g) 
(4.15h) 



The polynomials Cf for 9 < A; < 24 are reported in the Mathematica file data_FREE .m 
that is included with the preprint version of this article at arXiv.org; they can also be 
read off from the results of Section 5 below [cf. (5.25)/(5.33)]. Finally, the polynomials 
B^ for 9 < k < 24 can be determined from the C^ using (4.9). 

Actually, we can do better than this. We believe that the coefficients c^im) are, 
for each fixed k > 0, the restriction to integers m > m 
m of degree 1. If we compute the difference 



F 
mill 



(k) of a polynomial C^ in 



KM 



cl{m) 



CI (m) 



(4.16) 



between the numerical coefficients c|( 



m) and the corresponding polynomials C^, we 
find, not surprisingly, that they are nonzero whenever m < mjai^{k): see Table 3. 
If we could somehow guess an analytic form for at least some of these coefficients 
m), we could then deffne improved coefficients c^{Tn) by 



Af 



"fef 



cl[m] = c^ (m) — Aji. (m) , 



(4.17) 



so that these coefficients c^im) would be, for each fixed k, the restriction to inte- 
gers m > 'mJnjjj(A;) of the same polynomial C^, with a smaller threshold ^min(^) < 
^min(^)- The important point here is that a smaller threshold Tn^^:^^{k) implies that 
we can obtain more polynomials C^ with the same raw data. 

By inspecting Table 3, it is not difficult to realize that there are some patterns in 
A|(m) immediately below the threshold fn^:^j^{k): for instance, for odd k = 2p + 1 
with 3 < j9 < 10, we have AL_,_^(j9 — 2) = 1; and for even k = 2p with 4 < p < 10, 
we have A2p(j9 — 3) = 3p — 2. We then assume that this behavior holds true for all 



larger p. For other subsets of the nonzero values of A^ 



[m. 



slightly farther below the 



boundary m^-^j^{k), we likewise find simple polynomial Ansatze. Our empirical results 
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are: 



^l+lip- 


-2) = 


= 1 , for p > 3 


(4.18a) 


^l+lip- 


-3) = 


= f P^ - f p + 27 , for p > 4 


(4.18b) 


Aj,+i(p- 


-4) = 


= >^ 'lfp' + '\''p' 'Il'p + 917, forp>5 


(4.18c) 


KiP - 


-3) = 


= 3p-2, for p > 4 


(4.18d) 


K(p - 


-4) = 


= fp3-66p2+8iZp-l7, forp>5 


(4.18e) 


K(p - 


-5) = 


_ 473 5 2207 4 , 9023 3 38413 2 , 29763 , on/li 
" AOP 12 P ' 8 P 12 P ' W P^ ^^^^ ' 








for p > 7 


(4.18f) 



We are able to test each fit on at least one additional data point. Notice that in 
(4.18f), the condition p >7 does not follow the expected behavior from the previous 
correction terms (i.e., one would have expected p > Q). The new threshold ^min(A;) 
is given by 



"^mm(^) 



1 

2 

m 



5 



if < A; < 11 

if 12 < A; < 14 
if A; > 15 



(4.19) 



By this method, we can obtain the polynomials Cf (and therefore the polynomials 
B^) up to A; = 30. Indeed, the fits to obtain the polynomials Cf with A; < 30 were 
tested on at least one additional data point. 

If we do not demand to have at least one extra data point to test the fits, we 
can extend this computation up to A; = 32. We can then guess one further correction 
term A^(m) (this time with an additional data point to test the fit): 



^F f^ p;^ _ 161 



^2p- 



-i(p-5) = ^p' 



43693 5 
240 P 



21757 ^4 
16 ^ 



86909 ^3 



16 



po _^ 



110297 2 
8 i^ 



433772 
15 



p 



+42719 , for p > 8 . 
With this additional correction, the new threshold rnZ,: is 



(4.20) 



^mm(^) 



if < A; < 11 

if 12 < A; < 14 
if 15 < A; < 17 
if A; > 18 



(4.21) 



UlJ-5 

so the computation of the polynomials Cf can be extended up to A; = 33 (with no 
extra data points to test the fit). The polynomials Cf with 1 < A; < 33 are included 
in the mathematica file data_FREE.m. 
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4.3 Cylindrical boundary conditions 

We have likewise obtained the leading eigenvalue A^(m) for < -m < 13 as a 
power series in q~^ [cf. (4.3b)] through order k = 40. The resulting coefficients b^{m) 
are displayed in Table 4, and the corresponding coefficients c^{fn) are displayed in 
Table 5. As for free boundary conditions, we note that for all {k, m) that we have 
computed (i.e., 1 < m < 13 and < A; < 40), the coefficients b^{m) and kc^{m) are 
integers. 

We observe empirically that, for each fixed k, the coefficients b^{m) are the re- 
striction to integers m of a polynomial B^ in m of degree k, and that the coefficients 
c^{m) are the restriction to integers m of a polynomial C^ in m of degree 1, provided 
that we restrict attention to m> 'n^min(^) with 



m„ 



Xk) 



k + 2 . 



(4.22) 



Below this threshold, the coefficients deviate from polynomial behavior. The poly- 
nomial behavior can be extended downwards by two steps, i.e. to m = k, if we use 
the coefficients b^{m) [cf. (4.8)] in place of b^{m). With our available data together 
with the tricks described in the preceding subsection, we are able to determine these 
polynomials for < A; < 15. The coefficients b^{m) are displayed in Table 6. 

We begin, as before, by ffiting b^{m) for m > A; to a polynomial B^ of degree k. 
With our data we can do this for < A; < 6; we also have, in each case, at least one 
extra data point to test the ffi. Our results are: 



B^im 
B^{m 



5fl 



m 



m 



B 



[m 



BUm 



B!i 



m 



1 

2m 

2n? 



|m^ 



\m^ 



m 



2m2 
2m^ 



-,ra 



^rn^ + |m 



j^nv' - ^m^ + \m^ + ^m^ - ^m 



45'" 



|m^ + Im"^ 



^m3 



614^2 
45 '"' 



+ 



55 



m 



(4.23a) 
(4.23b) 
(4.23c) 
(4.23d) 
(4.23e) 
(4.23f) 
(4.23g) 



Notice that the three highest-order coefficients agree with those of the corresponding 
polynomial a^{m), i.e. 



B',{m) 



hlijn) 



for < A; < 3 



aXim) + Oim^ ^) for A; > 4 



(4.24) 



Note also that the constant term vanishes in all these polynomials except B^. 

As in the previous subsection, we can extract the desired information more ef- 
ficiently by analyzing the coefficients (^{m), which are found empirically to be, for 
each fixed A;, a polynomial C^ in m of degree 1 ior m > m^^^^{k) = k + 2, i.e. 



CUm) 



P nP 

aj.m + Pf. 



(4.25) 
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As we need only two coefficients for such a polynomial, we can obtain these polyno- 
mials up to /c = 9 (if we want at least one extra data point to test the fit) or A; = 10 
(if we don't). 

However, we can do slightly better if we consider instead of the coefficients c\ 



the modified coefficients c^{Tn) defined by (4.10). We find empirically that kc^{m) is 
an integer for all the computed values of {k, m): see Table 7. We also find empirically 
that c^{m) is, for each fixed k, the restriction to integers m > m^i^ik) of the same 
polynomial C^ with a smaller threshold m^^^{k) = max(A;, 2) < Tn^a:^^{k). In this way 
we can obtain the polynomials C^ up to /c = 11 or A; = 12, depending on whether or 
not we insist on having an extra data point to test the fit. Our results for k < 8 are: 



Cf(m 



C 



(m 






3(m 



a 



cl 



m 



C^{m 



,m 



C^{m 



(m 



-2m 



-m 



m 



m 



5 



b 



m 



m 



(4.26a) 
(4.26b) 
(4.26c) 
(4.26d) 
(4.26e) 
(4.26f) 
(4.26g) 
(4.26h) 



The polynomials Cl" for 9 < A; < 12 are reported in the Mathematica file data_CYL .m 
that is included with the preprint version of this article at arXiv.org; they can also be 
read off from the results of Section 5 below [cf. (5.25)]. Please note that the constant 
term vanishes in all these polynomials, while the term linear in m is the same as for 
free boundary conditions [cf. (4.15)]: 



a,. 



P! 



ai. 







(4.27a) 
(4.27b) 



in all cases that we are able to test (namely, 1 < A; < 12). Finally, the polynomials 
B^ for 7 < A; < 12 can be determined from the C^ using (4.10). 

These results can be improved in the same way as we did in the previous subsec- 
tion. First we compute the difference A^(m) 



cXim 



-cl^ 



m): 
p 

k 



see Table 8. We then 



try to guess an analytic form for some of the coefficients A^m), and we define im- 
proved coefficients c^(m) [as in (4.17)] so that the c^im) will be, for each fixed k, the 

^mm(^) of the polynomial C^ , with a smaller threshold 
^min(^) < "^min(^)- ^^ in the case of free boundary conditions, we find empirically 
that the coefficients A^{m) closest to the boundary ^nuii(A;) are the restriction to 



restriction to integers m > m^^^{k) of the polynomial C| 



35 



integers m of certain polynomials: 

Al{k-l) = i-lf-^i^k^-fk + A), for A; > 4 (4.28a) 

Al{k-2) = {-l)''-^fk^-10k^ + fk + l), for A; > 6 (4.28b) 

Al{k-3) = (-l)'=-3(fA;4-^P-i|lA;2-2||iA; + i9l) , for A; > 8 

(4.28c) 



p 

min 



Again, each fit can be tested on at least an extra data point. The new threshold mj^i 
is 

h iik<2 

^min(^) = UIJ + 1 if 3 < A; < 8 (4.29) 

A; - 3 if A; > 9 , 

By this method, we can obtain the polynomials C^" (and therefore the polynomials 
B^) up to A; = 14, with at least one extra data point to test the fit. 

If we do not insist on having an extra data point to test the fit, we can extend 
this computation of C^ up to A; = 15. We can then guess one further correction term 
A^{m) (again with no additional test for the fits): 

Al{k - 4) = (-1)'=-^ (f A;5 - i|5 A;4 + Z|lA;3 - IflA;^ - '-§k + 695) , for A; > 10 . 

(4.30) 
With this correction the new threshold is 

h iik<2 

^min(^) = UIJ + 1 if 3 < A; < 10 (4.31) 

A; - 4 if A; > 11 

so the computation of the polynomials C^" can be extended up to A; = 16 (with no 
extra data points to test the fit). 

The polynomials C^" for 10 < A; < 16 also have a zero constant term, and the 
relation (4.27) between periodic and free boundary conditions continues to hold. 

5 Thermodynamic limit (m — > oo) 

In previous sections we have dealt with semi-infinite square-lattice strips of fixed 
width m. In this section we will study the thermodynamic limit m ^ oo of the free 
energy of our model. 

In Section 5.1 we introduce some preliminary definitions and discuss the expected 
behavior of the strip free energies per unit length, /^(g) and f^{q), as a function 
of the strip width m. In Section 5.2 we discuss the large-|g| expansion of the bulk 
free energy. Bakaev and Kabanovich [12] have calculated the first 36 terms of this 
expansion; we confirm their computation and extend it by providing 11 more terms, 
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i.e. through order g~^^. We also compute the large-|g| expansion of the surface (resp. 
corner) free energy through order g~^'^ (resp. Q'"*^^). These computations are based on 
the finite-lattice method [44-49]. 

In Section 5.3 we obtain (using the polynomials Cf) a large-|g| expansion for 
the limiting free energy f^{q) of a semi-infinite strip of width m and free boundary 
conditions. We find that this expansion contains two terms: a bulk term (independent 
of m) and a surface term (linear in 1/m). This computation gives an independent 
check of the first 33 terms of the series expansions for the bulk and surface free energies 
for the square lattice. In Section 5.4 we repeat the computation with cylindrical 
boundary conditions. We find that the bulk contribution is the same as for free 
boundary conditions and that there is no surface contribution. These results provide 
a theoretical interpretation for some aspects of the behavior found in the preceding 
section for Cf and Cf [cf. (4.15)/(4.26) and (4.27)]. 

Finally, in Section 5.5 we perform a series-extrapolation analysis of the large-|g| 
series for the bulk, surface and corner free energies, in an effort to locate their singular 
points in the complex g-plane. 

5.1 Generalities and finite-size-scaling theory 

Corollary 3.3 shows that, for each width m and each boundary condition (free 
or cylindrical), the transfer matrix has, for sufficiently large |g|, a single dominant 
eigenvalue X^{q) that moreover is an analytic function of q (in fact, it is g™" times 
an analytic function of q^^). However, Corollary 3.3 gives no information about how 
large \q\ has to be for this behavior to occur; in particular, there is no guarantee that 
this large-g domain is uniform in m. However, the uniformity in m can be proven by 
invoking the following theorem: 

Theorem 5.1 [50, Corollary 5.3 and Proposition 5.4] Let G = (V, E) be a 

loopless^^ finite undirected graph of maximum degree A. Then all the zeros of the 
chromatic polynomial Pciq) ^^e in the disc \q\ < 7.963907A. 

Let us remark that this theorem has been recently improved by Fernandez and Pro- 
cacci [51, Corollary 2]: they showed that the constant 7.963907 in Theorem 5.1 can 
be replaced by 6.907652.^^ 

For the square lattice with any of the standard boundary conditions (free, cylin- 
drical, cyclic or toroidal) we have A = 4, so we can conclude that all chromatic roots 



^^ Warning: We are here using tlic grapli tlieorists' terminology, in which a loop is an edge 
connecting a vertex to itself. Obviously, if G has a loop, then Pciq) is identically zero. What 
physicists often call "loops" — particularly when referring to Feynman diagrams — are called 
"cycles" or "circuits" by graph theorists. Thus, a "3-loop Feynman diagram" is a graph with with 
cyclomatic number 3; it may or may not have loops. 

^^ Jackson, Procacci and Sokal [52] have recently observed that the Fernandez-Procacci constant 
is in fact 

K* = W{e/2)/[l-W{e/2)]^ w 6.907 651697 774 449 218... 

where W is the Lambert W function [53], i.e. the inverse function to a; h^ xe^ . 
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lie inside the disc |g| < 7.963907 x 4 = 31.855628. Actually, by [50, Corollary 5.3 and 
Table 1], we have for A = 4 the slightly stronger bound |g| < C(4) < 29.081607. With 
the improved result of Fernandez and Procacci, we get l^l < 6.907652x4 = 27.630607; 
and for the particular case A = 4 these authors [51, Corollary 1] obtained the slightly 
better bound |g| < C*(4) < 24.443218. 

It follows that the limiting curves Bm must also lie inside the disc |g| < 24.443218 
for all widths m: for if part of the limiting curve were to lie outside the (closed) disc, 
then the Beraha-Kahane-Weiss theorem (Theorem 2.2) would imply that chromatic 
roots would also lie outside the disc for mxn strips of all sufficiently large lengths n. 
Furthermore, using again the Beraha-Kahane-Weiss theorem we can conclude that, 
outside this disc, the transfer matrix for each width m must have one and only one 
eigenvalue of largest modulus. Since this dominant eigenvalue cannot collide with 
any other eigenvalue, it must be an analytic function of q outside the given disc. 

In summary, the transfer matrix for a square-lattice strip of width m and with 
free or cylindrical boundary conditions has a single dominant eigenvalue X^'m (q) that 
is an analytic function of q (in fact, g™ times an analytic function of q~^) whenever 
|g| > 24.443218. 

Let us now introduce the free energy per site for a finite strip with free or cylin- 
drical boundary conditions, 

and its limiting value for a semi-infinite strip, 

fTil) = lim— log PG.,(g). (5.2) 

Finally, let us introduce the free energy per site for the infinite lattice, 

/F/P(^) ^ li^ J- log P^(q). (5.3) 

m,n^oo rnn ^'^ 

Here we are assuming that the indicated limits exist and that in (5.3) the limit is 
independent of the way that m and n tend to infinity. Furthermore, it is natural to 
expect that in (5.3) the limiting free energy is independent of boundary conditions, 
in which case we can omit the superscripts F or P and write simply f{q). 

In fact, some of these assumptions can be proven. Indeed, the above discussion 
guarantees that at least for |g| > 24.443218, the limiting strip free energy fm{q) exists 
for all m and is given by 

fTiq) = -^ogXZiq), (5.4) 

m 

which in particular is an analytic function of q in the indicated domain. Moreover, 
Procacci et al. [54, Theorem 2] have proven that, when |g| is large enough (namely, 
|g| > 8e'^ ~ 160.684295), the infinite-volume limiting free energy f{q) exists and 
is analytic in 1/g and is the same for all reasonable sequences of graphs Gmxn (in 
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particular, it is the same for free, cylindrical, cyclic and toroidal boundary conditions 
and is independent of the way that m and n tend to infinity). In this paper we will 
take n ^ oo first and then take m — *> oo, so that 

f{q) = hm /^/P(g) . (5.5) 

m— >oo 

Finite-size-scaling theory [55, Section 2.5] gives a rather precise prediction for the 
form of the free energy (5.1)/(5.2) for a finite or semi-infinite system away from a 
critical point (and in the absence of soft modes). In particular, for an m x n strip 
with free or cylindrical boundary conditions and bulk correlation length ^buik *^ ^; n, 
the predicted behavior is 

fin = /bulk + ^^^^/surf + — /corner + 0(e-""°("^'")/«--) (5.6a) 

mn mn 

fZ,n = /bulk + -/surf + 0(e-""°('"'")/«--) (5.6b) 

n 
where /bulk = /, /surf and /corner are, respectively, the bulk, surface and corner free 
energies. (More precisely, /surf is the free energy for two units of surface, and /comer 
is the free energy for four corners.) For a semi- infinite strip m x oo with free or 
cylindrical boundary conditions, we have 

fl = /bulk + ^/surf + 0(e"'"/«'^"'>^) (5.7a) 

fl = /bulk + 0(e-"^/«--) (5.7b) 

The relations (5.7) of course hold for the chromatic polynomials at fixed large q. 
But we can also argue heuristically what they should imply for the series expansion 
in powers of 1/g. It is not difficult to see that, for large q, we have 

e-i/Cbu,kW = I ^ o(\) (5.8) 

q \q^ / 

(just as for a one- dimensional Potts antiferromagnet at zero temperature).^'^ We can 
therefore interpret 0(e~™/^''""'*^''^) as meaning 0{q~'^). Therefore, we expect that 

fl = /bulk + -/surf + 0(g-™) (5.9a) 

m 

fl = /buik + 0(g-™) (5.9b) 



^°Let G — (V, E) be a finite graph (let us suppose for simplicity that it has no loops or mul- 
tiple edges) with n vertices and m edges; then one sees immediately from the Fortuin-Kasteleyn 
representation (2.4) that its Potts- model partition function has the large-g expansion 

Zaiq.v) = q- + mq-' + ^^^^ ^) ^^-2 ^ ^(q"-^) 

When G is a finite piece of a regular lattice, the corresponding expansion for |T^|~^ log Zciq, v) gives 
in the infinite-volume limit the lavge-q expansion for the bulk free energy f{q,v). 

Now let i,j be vertices of G; then the unnormalized 2-point correlation function ZaicTi ■ ctj) = 
ZQ{{q5ij^^„, — \-)/{q — 1)) is given by a representation like (2.4) but with the constraint that i and j 
must belong to the same connected component. The dominant terms of the large-g expansion are 
the ones in which this component has the minimal number of vertices and all other components are 
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This explains why the coefficients c^ (m) in the expansion of logA^^m, = rnfm 
[cf. (4.9)] are polynomials of degree 1 in m for large enough m, i.e. m > rnj^^{k), 
and consequently why the coefficients 6^ (m) in the expansion of Ay,,m [cf. (4.3)] 
are polynomials of degree k in m. It also explains why, for cylindrical boundary 
conditions, c^{m) is strictly proportional to m, i.e. the constant term vanishes [cf. 

P /P 

(4.26)/(4.27b)]. The error term 0(g^™) in (5.9) furthermore suggests that m^^^ [k) = 
k+1. Of course, we should not take too seriously the "+1" here, since the amplitude of 
the correction term in (5.7) could be proportional to a positive or negative power of q. 

P /P 

But we do predict that m^^^^ (k) = k + 0(1) as A; — > oo. This is indeed what we found 
for cylindrical boundary conditions, where we have iTL^^^^{k) = k + 2 [cf. (4.22)]. For 
free boundary conditions, however, we found the faster convergence rn^^^^{k) ^ k/2 
[cf. (4.11)], for which we lack at present any theoretical explanation. 

In the next subsection we will use our transfer matrices to compute the large-g 
expansion for the bulk, surface and corner free energies, using the finite-lattice method 
[44-49]. In the subsequent two subsections we will compute the large-g expansion of 
the strip free energy fm{q) for free and cylindrical boundary conditions, respectively. 



isolated vertices; one therefore gets 



^G^O-i ■ CTj 



/ ?; \ "G V^-j ) 

Qg{v;i,j)[-) g"[l + 0(g-i)] 



where dQ{i,j) is the length of the shortest path in G from i to j, and Qg{v] i,j) is a polynomial in 
V that enumerates the connected subgraphs of G that contain i and j and have exactly dQ{i,j) + 1 
vertices [with a weight v for each edge beyond the minimum number dQ{i,jy\. In particular, if G is 
triangle- free, these subgraphs are simply shortest paths from i to j. In general Qg{v] i,j) can grow 
exponentially in dQ{i,j) [when e.g. G is an infinite regular lattice]; but if G is a piece of the square 
lattice and i — j lies along an axis direction, then QG{v;i, j) = 1. It follows that the exponential 
decay rate of correlations along an axis is 



g-i/Cbuik(g,f) _ 



fJ- 

VloP 



o 

' 19 



For V = —1 this gives (5.8). 

It is instructive to ask how these results would be seen in the transfer-matrix formalism. Ordinarily 
one has g-i/ibuik — |A2/A^|, where A^, is the dominant eigenvalue and A2 is the first subleading 
eigenvalue. But when one performs the computation using our transfer matrices, one finds A2/A^ = 
am {v)/q'^ + 0{l/q'^) for suitable polynomials am — not the predicted v/q. (Indeed, for width 
m = 1 — i.e., a one-dimensional Potts model — the transfer matrix is of size 1x1, i.e. there is no 
subleading eigenvalue at all.) What is going on here? 

The point is that the exponential decay rate in the correlation function {ctq ■ cTx) is controlled 
by a "colored" intermediate state, i.e. the state obtained by applying the field erg to the vacuum. 
The corresponding eigenvalue A2 would be seen in a transfer matrix in the spin representation [6, 
Section 3.1]; but it is not seen in our transfer matrix in the Fortuin-Kasteleyn representation, which 
represents only "colorless" states (i.e., states invariant under the Potts global symmetry group Sq). 
Rather, the first subleading eigenvalue of the latter transfer matrix corresponds to a "colorless" 
two-particle state, hence has A2/A* of order 1/g^. More precisely, for square-lattice strips of widths 
771 > 4 (resp. TO > 3) with cylindrical (resp. free) boundary conditions, we find (at least up to 
TO, = 9p, 7f) that there are at least two subleading eigenvalues of order l/q'^. Some of these satisfy 
A2/A* = {v/q)-^ + 0{q~^), and the others satisfy A2/A* = {1 + v){v / q)-^ + 0{q~^); note that the latter 
ones vanish at order 1/q'^ for the chromatic polynomial v = —1. All other eigenvalues are 0(q^^). 

40 



using our polynomials Cf and C^ . 



5.2 Large-g expansion for the bulk, surface and corner free 
energies 

First of all, as we are interested in the large-g limit, it is convenient to explicitly 
remove the leading term log q in the free energy by considering the modified chromatic 
polynomial Pq for a loopless graph G = {V, E): 



Using the Fortuin-Kasteleyn representation (2.4) we get 

Pciq) 



\^(_l)l^lgfc(^)-|V| 



ACS 



E(-i)""(v?) 



\A\~c{A) 



(5.10) 

(5.11a) 
(5.11b) 



ACE 



where c{A) = \A\ — \V\ + k{A) is the cyclomatic number of the subgraph (V, A). 

It is instructive to begin by computing "by hand" the first few terms of the 
large-g expansion for the bulk, surface and corner free energies. To do this, let us first 
consider an m x n square lattice with free boundary conditions: it has \V\= mn sites, 
\E\ = 2mn—m—n edges, and |F| = (m— l)(n— 1) square faces ( "plaquettes" ) . We can 
compute the first few first terms in the large-g expansion for the modified chromatic 
polynomial Pmpxnpiq) by using (5.11b) and explicitly identifying the subsets A having 
a given small value of |A| — c(A): 



Only A = 0. 
A = any single edge. 
A = two distinct edges. 
A = three distinct edges or four edges forming a 



\A\-ciA) = 

\A\-c{A) = l 

\A\-ciA) = 2 

\A\ - c{A) = ; 
plaquette. 

\A\ — c{A) = 4: A = four distinct edges not forming a plaquette or 
four edges forming a plaquette together with one additional edge. 



We therefore have 



^mpxnp [q) 



1- 



E\ \E\(\E\- 
2g2 



q 



\E\{\E\-l){\E\-2) 



6 



+ 



\E\{\E\-l){\E\-2){\E\-3) 



24 



-\F\i\E\-3) 



+0{q-') . 



(5.12) 
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Taking the logarithm, dividing by \V\, and putting back the leading term logg, one 
finds the large-g expansion for the free energy 



J m,n\^/ 



2 115 

logg ^ + TT + TT 

q g^ Sq-^ 2q^ 



1 



2 
3^ 



11 



1 1 
— + - 
m n 



+ 



ri 31 


1 


y + ^. 


mn 



0{q 



-5n 



Comparing to the finite-size-scaling Ansatz (5.6), we obtain 



/buik(g) 

/surf(g) 
/corner W J 



2 1 1 

logg ^ + TT 

q g^ 3g'^ 



2? + «(«"' 



1 

- + 

1 



1 

V 

3 

g3 g4 



2 

3^ 
0(g 



11 



0{q 



-5\ 



-5\ 



(5.13) 

(5.14a) 
(5.14b) 
(5.14c) 



Now consider an m x n square lattice with cylindrical boundary conditions: we 
have \V\ = mn, \E\ = 2mn — m, and \F\ = nfi{n — 1). For small m we have to 
worry about terms A that wind horizontally around the lattice using the periodic 
boundary conditions; such terms start at order g~(™'^i). But we can avoid such terms 
simply by assuming that m is large enough, i.e. m > A; + 2 if we want an expansion 
valid through order q~^. Therefore, we can obtain the expansion through order g~^ 
by assuming that ?7i > 6; then the contributing terms are exactly the same as those 
for free boundary conditions, and we obtain the same expansion (5.12) but with the 
modified values of \E\ and \F\. A simple computation shows that 



fm,n{^) = /bulk(g) + /surf(g) " + ^(g" 

n 



(5.15) 



where /bulk and /surf are the same as those given in (5.14), in agreement with the 
finite-size-scaling prediction (5.6). 

Unfortunately, this elementary graphical method for generating the large-g expan- 
sion is not very efficient if one wants to go to high order. To obtain long expansions we 
will use a more sophisticated procedure: the finite-lattice method pioneered by Enting 
and collaborators [44-49]. Indeed, Bakaev and Kabanovich [12] used this method to 
obtain, already in 1994, the large-g expansion of the bulk free energy through order 
g,-36 21 'jj^gij; results can be summarized as follows: as |g| -^ oo, the exponential 
of the bulk free energy per site for the zero-temperature Potts antiferromagnet (i.e., 
chromatic polynomial) on the square lattice is given by the series expansion 



J{q) 



{<! 



'-[l + z^ + z'^ + 3z^ + Az^ + 3z^° + 3z" + llz^^ + 24z^^ + 8z^^ 
-91z^^ - 261^^^ - 290^^^ + . . . - 3068121066^^ + 0(/^)] (5.16) 



2ipQ]. relevant earlier work on obtaining large-g series for the infinite- volume limit of the chromatic 
polynomial, see Kim and Enting [46] and the references cited therein. 
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(the full series is given in Table 9), where z is defined as 

1 



g-1 



(5.17) 



Before presenting our extensions of this result, let us first briefiy review the finite- 
lattice method. 

In the finite-lattice method [44-49], the large-g expansion of the infinite- volume 
free energy through a given order in z = l/(g — 1) is written as a linear combination 
of free energies for rectangles of various sizes r ^ s: 

/(g)-logg = lim logP^^y^nA^) (5.18a) 

m,n-»oo mn 

= Yl Mr,s)logPr,^s,{z) + 0{z'''-'), (5.18b) 

{r;s)eB{k) 

where we have used the variable z = l/(g — 1) instead of g = 1 + 1/ z, and the modified 
partition function PmpxriF [see (5.10)]. The sum in (5.18) is taken over all rectangles 
r X s belonging to the set 

B{k) = {{r,s): r < s and r + s < k} . (5.19) 

The weights ak{r, s) are defined as 

J2Wk{r,s) forr<s 
Mr,s) = <,^, , , , (5.20) 

I Wk{r,r) tor r = s 

where 

1 for r + s = k 

—3 for r + s = k — 1 

Wk(r,s) = <3 forr + s = A;-2 (5.21) 

—1 for r + s = k — 3 

otherwise 

The error term in (5.18) is given by a particular subclass of connected graphs that do 
not fit into any of the rectangles in B{k) [44] [56, Chapter 12]. In our case, these graphs 
are [12] the convex polygons^^ of perimeter 2/c — 2 having the property that any pair of 
nearest-neighbor sites belonging to the polygon must be connected by an edge of the 
polygon. ^^ Since any polygon of perimeter 2A; — 2 has |A|— c(A) = (2A;— 2) — 1 = 2A; — 3, 
we deduce that the error term in (5.18) is of order g~(2fc-3) ^ ^2fc-3^ 

Now suppose that we are somehow able to compute the transfer matrices for widths 
L < Ljnax with free boundary conditions. We can then use these transfer matrices to 



22\Ye recall [57] that a convex polygon (in the square lattice) is a self-avoiding polygon whose 
length equals the perimeter of its minimal bounding rectangle. 

^•^This condition excludes, for instance, polygons containing a "bottleneck" of width 1 or a "pro- 
tuberance" of width 1. 
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compute the partition functions Prpxsp for r < Lmax and s arbitrary. This means that 
if we set the cut-off k equal to 2Lmax + 1, we will be able to compute the partition 
function for all pairs (r, s) G B{k). It follows that formula (5.18) gives the bulk-free- 
energy series correct through order ^4Lnmx-2^ -^^ have empirically checked (by doing 
computations for different values of Lmax and checking to what order they agree) that 
this formula is correct. In our case Lmax = 12, so we expect to obtain the free-energy 
series correct up to order z"^^. 

If we compare the results coming from different values of -^max, we find empirically 
that the first incorrect term (of order ^^Lmax-i^ jg given by the generating function of 
the aforementioned subclass of convex polygons, which is [12] 

>^ 2k _ 8 2-2a:'-a:Vl-4x2 ,^ 3 - Ax^ - AVl - Ax^ 

^/''^ " "^ (l-4a;2)(2 + a;2) + "^ {1 - Ax^^ ' ^^'^^^ 

If we add the term /i2Lmax-2^'^'"'"'~^ ^o ^^e series obtained with a given value of L^a.x, 
we get a series correct up to (and including) ^4-^max-i_ j^ q^j, ^g^gg 2,^^^ = 12, so we 
obtain a series expansion correct up to order 2;^^. The idea of adding this correction 
term is due to Bakaev and Kabanovich [12] and was used by them to improve their 
own series by one term.^^ 

In this way, we confirm the expansion (5.16) found by Bakaev and Kabanovich [12] 
and extend it by 11 terms: 

g/W ^ liull! [1 + ^3 ^^7 ^3^8 ^4^9 ^3^10 ^3^11 ^^^^12 ^24^13 ^8^14 
q 

-91z^^ - 261z^^ - 290^^^ + . . . -598931311074/^ + 0(/^)] , (5.23) 

where the full series is given in Table 9. In terms of the variable 1/g, we obtain 

g/(g) = g [1 _ 2q-^ + g"2 + q'^ + g"^ + q'^ + q'^ + 2q-^ + 9g"^ + 38g'^ + 130g"^° 
+378g-" + 987g-^2 + . . . + 1311159363081366872g-^^ + 0{q-^^)] . (5.24) 

Finally, for the bulk free energy /(g) itself (rather than its exponential) in terms of 
the variable 1/g, we obtain 

1 I 1 I 5 I 28 I 55 I 89 I 81 I 505 I 1029 , 7742 



f(o) - loe-rJ-^-i + ^- + -2- + ^2._L.55_,J9_,^^, 505, 1029, 

25291 _|_ 73552 , 197755 , , 19001827661! 
12gi2 "T I3gl3 "T ^4gi4 -|- • • • -|- 47g. 



q q'^ ' 3i3=* ' 2ij-* ' bq^ ' 6qf* ' Jq'' ' 4q^ ' 9q^ ' 5g"' ' llg" 
, 25291 I 73552 , 197755 , , 190018276619486037135 , nfr,-'iS\ (K 0K\ 



Clearly, e^^'^^ has a much simpler expansion than /(g); in particular, its coefficients 
are integers (at least through the order calculated thus far).^^ A further simplification 
is obtained by using the variable z = l/(g — 1) in place of 1/g: the integer coefficients 



^"^They actually went farther and improved their series by a second term, by enumerating a more 
complicated class of contributing graphs that we refrain from considering here. 

^^The coefficients of f{q) are not integers, but k [q^'']f{q) is an integer (at least through the 
order calculated thus far). Indeed, it is not hard to show that if F{z) is a power series with integer 
coefficients and constant term 1, then k [z''] logF(z) is always an integer. 
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become much smaller. Finally, a slight extra simplification arises from extracting the 
prefactor (g — l)^/g in (5.23). 

The finite-lattice method can be used for computing series expansions, not only 
of the bulk free energy, but also of the surface and corner free energies [cf. (5.6)]. The 
main ideas were introduced by Enting [45] three decades ago and have been further 
developed by him [48,49]. To our knowledge, however, no one has yet applied this 
technique to the large-g expansions of chromatic polynomials, so the all the results 
to be presented below are new. 

As for the bulk case, the large-g expansions of the surface and corner free energies 
through a given order in z = l/(g — 1) are written as certain linear combinations of 
free energies for rectangles of various sizes r x s [cf. (5.18)]: 



/surf(g) = 5Z Pk{r,s)\0gPrpy,sp{z)+0iz 

{r,s)(^B{k) 

/comer(g) = JZ 1 k{r, s) log Pr^xspiz) + 0{z 

(r,s)G_B(fc) 



2k-3\ 



2k-3\ 



(5.26a) 
(5.26b) 



where the set B{k) is the same as for the bulk case [cf. (5.19)]. The weights for the 
surface free energy are given by [45] 



Pk{r,s) 



Sk{r,s) + Skis,r) 
Sk{r,r) 



for r < s 
for r = s 



(5.27) 



where 
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The weights for the corner free energy are given by [45] 

lk{r,s) 
where 

VJr-s) = 



. \2Vk{r,s) 


for r < s 






'^ [Vk{r,r) 


for r = s 
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otherwise 







(5.29) 



(5.30) 



The error term in (5.26a,b) is given by the same reasoning as for the bulk case, i.e. 
it is 0{z'^^~^). Therefore, by considering strips up to a maximum width Lmax, we 
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can use (5.26a,b) to obtain the series for the surface and corner free energies correct 
through order z'^^'^'^^"'^ . We have empirically checked (by doing computations for 
different values of -^max and checking to what order they agree) that this formula is 
correct. In our case we have Lmax = 12, hence series valid through order 2;^^. 

For the surface free energy we can conjecturally extend the series by one term 
by comparing the series for different values of Lmax and noticing empirically that the 
needed correction for the term of order ^^Lmax-i jg simply a multiple of the correc- 
tion term (5.22) for the bulk-free-energy series. More precisely, if we add the term 
— -^maxAi2Lmax'2^'^'"'''' ""^ ^^ the serics obtained with a given value of ivmax, we get a series 
correct through order ^^Lmax-i^ j,^ q^j~ ^g^gg L^^x = 12, so we obtain a series expan- 
sion for /surf correct (conjecturally) through order 2;^^. Unfortunately, we have not 
yet succeeded in figuring out an analogous correction for the corner free energy. 

The results for the series expansions are given by 



p/surf _ 


= 1 + z-z^-z^ + z^ -z'^ -8z^ - 16z^ - IQz^^ - I2z^^ - Alz^^ 




-138z^^ + ...- 130312353695974^^^ + 0{z^^) (5.31) 


^/corner 


= l + z^ + Az'^ + 12/ + 20/ + 28^1° + 67^1^ + 208^^^ + ASAz^^ 




+753z^* + ... + 448320847685638/^ + 0(/^) , (5.32) 



where the full list of coefficients is displayed in Table 9. In terms of the variable 1/q 

we get 

f _ lij^ 2 n 29 28 ^_243_ 1019 _ 4489 _ 17280 

^surf g "•" 2^2 3g3 4g4 5g5 3^6 7^7 g^s 9^9 lOgiO iigll 

14654 183143 , 1103009229135728011786 , n('^-48\ /cr oo\ 

3gi2 13<?13 + • • • 47(?47 + <^ly ) [0.66) 

f — J__I_A_I_A_I_J^_I_16_|_101_|_685_|_948_|_ 3409 

J comer — ^3-^54-^^5-^ 2g6 ^ q' ^ 2q» ^ 3q^ "^ gl" "^ g" 

+ f^ + ^ + . . . + 280315319437238591517 + ^(g-^^) (5.34) 

5.3 Large-g expansion for the strip free energy: 
Free boundary conditions 

We can check the results (5.25)/(5.33) for the bulk and surface free energies by 
computing the large- |g| expansion of the limiting free energy for a semi-infinite strip 
with free boundary conditions: 

1 
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\ogq + -y {-lfcl{m)q-' (5.35) 



m ^—^ 

fc=i 



If m > 'n^min(/i;), we can replace the coefficients cf (m) by the corresponding polyno- 
mials C^{m) [cf. (4.15)]. These polynomials Cf are of degree 1 in m for 1 < A; < 33, 
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and we have conjectured that this behavior holds for all values of A; > 1. Thus, fmio) 
contains [as predicted in (5.9a)] only two terms: a bulk term /buik(?) = fXl) that is 
independent of m, and a surface free energy energy ff^^^^q) that is of order 1/m: 

fl{q) = /buik(g) + -/surf(g) • (5.36) 

m 

Here = denotes that the two sides agree at each order q^^ of the expansion in pow- 
ers of q^^, but we require this only for m > iTLf^:^j^{k) . This computation using the 
polynomials C^ thus provides an independent check of the first 33 terms of the series 
(5.25) and (5.33). Conversely, the finite-lattice computation of the preceding subsec- 
tion provides an independent confirmation of the empirically observed regularities in 
the behavior of the polynomials Cf [cf. (4.16)-(4.21)]. 

We can slightly extend the check for the surface-free-energy series by using the 
series (5.25) for /bulk as an input: in this case, each polynomial Cf contains a single 
unknown coefficient to be determined (rather than two unknown coefficients). We 
then obtain the coefficient of the term z^"^ in /surf, which agrees with the result from 
the finite-lattice method displayed in Table 9. 

5.4 Large-g expansion for the strip free energy: 
Cylindrical boundary conditions 

The large- 1 g| expansion of the free energy for a semi-infinite strip with cylindrical 
boundary conditions is 



/:(?) = 


= -logAr,„(g) 
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= logg + — log 
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)_^{-l)%{m)q->^ 

.k=0 
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= logg + -V( 
k=i 


-lM{m)q-' 



(5.37) 

Once again, if m is large enough (depending on A;), we can replace the coefficients 
c^{m) by the corresponding polynomials C^{m) [cf. (4.26)]. These polynomials C^ 
are of degree 1 in m with a zero constant term for 1 < A; < 16, and we conjecture 
that this behavior holds for every A; > 1. Thus, /„(?) contains [as predicted in (5.9b)] 
only a single m-independent term fhuitiq) = f{q)] there is no surface contribution. 
Of course this result is to be expected, as an infinitely long cylinder has no boundary. 
This computation using the polynomials C^ yields a series for /bulk that is the same, 
up to the order we are able to compute (namely, q~^^), as for free boundary conditions; 
in particular, it provides an independent check of the first 16 terms of the series (5.25). 
Conversely, the finite-lattice computation of the preceding subsection provides an 
independent confirmation of the empirically observed regularities in the behavior of 
the polynomials C? [cf. (4.28)-(4.31)]. 

If we use the bulk free-energy series (5.25) as input, we can compute the term of 
order g~^^ in the surface free energy: as expected, we find that it vanishes. 
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5.5 Analysis of the large-g series 

In this subsection we shall perform a series-extrapolation analysis of the large-g 
series (5.23) ff. for the bulk free energy, (5.31) ff. for the surface free energy, and (5.32) 
ff. for the corner free energy. As a warm-up, we shall first perform an analogous 
analysis for the large-g series of the strip free energy fmig) for two selected cases, 
m = 3f and m = 4p. In these cases we can easily generate long large-g series (up to 
100 terms or more) using the symbolic power method (see footnote 16 above) and can 
compare the predicted location and nature of singularities with the exactly known 
answers [6]. In this way, we can learn how many terms in the large-g series are likely 
to be needed in order to extract specific features of the singularity structure. 

Before beginning this analysis, it is useful to know what types of singularities we 
are expecting. For statistical-mechanical models with only one dimension tending to 
infinity (i.e., finite-width strips), the answer is clear: the free energy fmig) [cf. (5.4)] 
is an analytic function of q except at branch points arising from the collision of two 
eigenvalues. Therefore, the limiting curve Bm is not a curve of singularities of the 
free energy; only its endpoints are singularities. 

For statistical-mechanical models with two or more dimensions tending to infinity, 
as in our infinite-width limit m —^ oo, the situation is much less clear. It is well 
established in a variety of cases [58-61] that the infinite-volume free energy has a 
"soft" essential singularity (i.e., one in which the free energy is infinitely differentiable 
but not analytic) at a first-order phase transition in the physical region. What is less 
clear, however, is how the free energy behaves on complex phase boundaries (such 
as our limiting curve Boo). For instance, the point h = is known rigorously to be 
an essential singularity for the free energy of a low-temperature Ising ferromagnet 
in dimension d > 2 [58,59,61]: it is not possible to analytically continue through 
this point, because the derivatives at /i = grow too fast (the Taylor series has zero 
radius of convergence). But might it be possible to analytically continue around this 
point in the complex /i-plane? If so, then the imaginary axis Re/i = 0, Irak ^ 
would be merely the curve where two analytic functions defined in a neighborhood 
of that curve — namely, the free energies that are initially defined for Re /i > and 
for Keh < but are analytically continuable to at least part of the opposite half- 
plane, the point h = excluded — have equal real part. By contrast, the alternative 
scenario is that these two analytic functions cannot be continued into any part of the 
opposite half-plane: that is, they would have the full imaginary axis Re /i = as a 
natural boundary. It does not seem to be known, even heuristically, which of these 
two scenarios is correct. ^^ 

5.5.1 Differential approximants 

To perform the series-extrapolation analysis of our large-g series, we shall use the 
differential- approximant (DA) method: see the review [64] and the references cited 
therein. The idea is simple: The i^th-order approximant [iVg, A^i, . . . , A^i^-; M] to 
a power series F{z) = Yl!k=ofkZ^ ^^ built by choosing polynomials Qo.Qi^ . . . ,Qk 



^^See also [62,63] for related work on complex phase boundaries and partition-functfon zeros. 
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and P of degrees No,Ni, . . . , Nk and M, respectively, so that the solution F of the 
inhomogenenous linear differential equation 



EQ.w(^^)'i^W = P{^) (5.38) 



j=0 

agrees with the initial coefficients of the series F{z) through order z^ , where A^ = 
X]j=o(-^* + 1) +M (notice that we fix the overall normalization with Qk{0) = !)• The 
singularities of F{z) are located at the zeros {ze} of the polynomial Qk (and possibly 
also at 2; = and z = 00), and the critical exponent associated to a simple zero z^ of 
Qk is given by 

A, = ir-l-%iM, (5.39) 

in the sense that F(z) ~ (2; — -2^)''*'^ as z ^ Z£. 

In practice, we have used a Mathematica code to obtain the polynomials Qj 
exactly (i.e., with exact rational arithmetic) from the coefficients {fk}, and then used 
the program MPSOLVE [65,66] to compute the N^ zeros of the polynomial Qk to 
very high precision (i.e., 100 digits). Our code and the numerical accuracy of its 
results were checked in a previous work [67] with the help of an independent C-\ — h 
program written by Y. Chan, A.J. Guttmann and A. Rechnitzer. 

In many cases, one is principally interested in the singularities that are located on 
the real axis. This assumption simplifies a bit the practical procedure, as one needs 
to search only on a one-dimensional space. However, in our case we want to locate 
all the singularities (real or complex) of the approximants. In fact, we expect from 
the sign pattern of the series coefficients (5.23)-(5.25) [cf. Table 9] that the leading 
singularity of the bulk free energy is not on the real z-axis. Therefore, we need to 
slightly change the usual definition of a non-defective approximant (i.e., those that 
are taken into account in the computation) [64,67]. Roughly speaking, we want to 
consider an approximant [Nq, Ni, . . . , Nk', M] to be non-defective if there is a zero of 
Qk sufficiently near to the expected (complex) singularity, and this zero is sufficiently 
well separated from all other zeros of Qk- 

Of course, this definition needs to be made more precise; and things are made more 
difficult by the fact that we do not know even roughly the position of the singularities 
in the complex z-plane. In order to find a rough estimate of these positions, we first 
make a histogram of all the "well-separated" complex zeros {ze} coming from all 
computed approximants [A^o, A^i, • • • , Nk] M]. A pair of complex zeros zi, z^ of Qk 
is defined to be well-separated \i\z\ — Z2\> R for some (small) free parameter i?; a 
zero Zx is then defined to be well-separated if it is well-separated from all other zeros 
(in particular, it is required to be simple). Given a fixed value of i?, we select from 
each polynomial Qk the well-separated zeros and make the corresponding histogram, 
using square cells of side B! = 0.04. This histogram is expected to display peaks at 
the singularities of the function F. We further expect that the final result should 
not depend strongly on the chosen value of R. In practice, we started with R = 0.2 
and then compared the histogram to those obtained with smaller values oi R = 
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0.1, 0.05, 0.02, etc. We considered that the procedure had converged when the number 
and positions of the peaks did not vary when "halving" the value of R. Usually, the 
choice R = 0.05 was optimal. 

Each peak displayed in the histogram corresponds to a cell in the complex z- 
plane: Re^ G [^1,^2] and Im^ G [Bi, B2] with A2 - Ai = B2 - Bi = R' = 0.04. 
Therefore, for a given singularity ze, we consider an approximant [Nq, Ni, . . . , Nk] M] 
to be non-defective if there is a (simple) zero of Qk inside the above cell, and there 
is no other zero of Qk inside the larger region Rez G [Ai — 0.05,^2 + 0.05] and 
Iraz G [Bi — 0.05, i?2 + 0.05].^^ Then we perform the analysis explained by Guttmann 
in his review article [64] . Taking all the non-defective approximants of a given order 
A^ corresponding to a given singularity Z£, we form the mean and standard deviation 
for each of the three quantities Kez, Imz and A. Of course, this standard deviation is 
only pseudo-statistical; for brevity we term it the "error bar" , but it is not necessarily 
indicative of the accuracy of the estimate. 

Our initial goals are twofold: to test that the above-described procedure gives the 
right answers in exactly soluble cases, and to determine how the estimates for the 
positions {zi} and critical exponents {A^} of the singularities vary with the available 
order N of the series. ^^ In order to achieve these goals, we have chosen two simple 
examples for which we know the position and type of all the singularities, and for 
which we can easily obtain many terms of the corresponding free energies (say, 100 
terms): namely, the strip free energies m = 3f and m = 4p. As in Ref. [67], we 
have computed the approximants of first and second order {K = 1,2) satisfying 
\Ni — Nj\ < 1. There are four free parameters in this computation: the maximum 
degree A^o.max of the polynomial Qq, the maximum degree M^ax of the inhomogeneous 
term P, the minimum number of series coefficients used in the analysis A'^rniru and the 
number of available coefficients N. We have chosen the first three parameters in terms 
of N in the following way: A^o.max = 0.8A^, Mmax = 0.5A^, and A^min = 0.3A^. Thus, 
the only free parameter in the computation is A^. In our test cases, we have chosen 
A^ to range between 30 and 100 in steps of 10. The two test cases are described in 
the next subsection. 

5.5.2 Test cases: m = 3f and m, = 4p 

The eigenvalues for the case m = 3f are the solutions of the quadratic equation [6] 

A^ - (g^ - 5g^ + llg - 10) A + (g^ - 7q^ + 19g^ - 24g + 11) = . (5.40) 

Therefore, the only singularities are those corresponding to the six endpoints of the 
limiting curve B3, where the discriminant of the quadratic vanishes. In the complex 



27pQ]. j.gg^i singularities, we used a slightly different definition: namely, for zeros near z = 0.5 (resp. 
z — —1) we asked that there should be no other zero in the region | Imz| < 0.1 and Rez G [0, 1.2] 
(resp. Reze [-1.8,0]). 

^^That is, we assume that we know the coefficients of the series F{z) through order z^ . Sometimes 
we shall refer loosely to N as the "number of available coefficients" even though, strictly speaking, 
the number of available coefficients (including the constant term) is N + 1. 
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plane oi z = l/{q — 1), they are located at 

zi ^ -0.2811172691 ± 0.7752009092 i (5.41a) 

Z2 ^ 0.4477839366 ± 0.5382490441 i (5.41b) 

Z2 = ^ + ^^^ ~ 0.5 ± 0.2886751346 i (5.41c) 
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In Figures l(a,b) we plot the histograms obtained with A^ = 50 and A^ = 100, 
respectively, using R = 0.05. The number of zeros contained in each square cell of 
linear size 0.04 is indicated with a gray scale of 16 tones: those cells with the maximum 
number of zeros are depicted in black, while those with the smallest number of zeros 
(i.e., less than 7% of that maximum number) are depicted in white. We have labelled 
the peaks with Im^; > as zi, . . . ,zj. We observe several empirical properties: as 
A^ is increased, the number of peaks increases (e.g. zq and z-i appear for N = 100 
but not for A^ = 50) and some of the peaks become sharper (e.g. Z4); furthermore, 
we find more peaks than the actual number of singularities (we expect three pairs of 
complex-conjugate peaks zi, Z2, Z3, but we also find two further pairs Z4, zj of complex- 
conjugate peaks and two peaks Z5, zq on the real axis).^^ For each of the seven peaks, 
we have performed the analysis described in the preceding subsection, in order to 
estimate the position of the singularity and the corresponding exponent. In Tables 10 
and 11 we display the results. 

We indeed find the three "correct" singularities associated to each of three end- 
points of the limiting curve in the upper half-plane (see Table 10). As the number of 
coefficients is increased, the accuracy of the results increase: with A^ = 100 we attain 
a maximum accuracy of order 10~^ (resp. 10~^) for the position of the singularity 
(resp. the critical exponent), and in all cases the estimate agrees with the known 
exact value within 3 times the claimed error bar (and usually much less). 

In Table 11, we show the "singularities" found in the series analysis that do not 
correspond to any of the endpoints of the limiting curve. Once again, as we increase 
the number of coefficients A^ in the analysis, the claimed error bars decrease. Let 
us stress that the exponents of the first three singularities (^4,^5,2:6) are consistent 
with A = 0; therefore, the DA analysis is asserting (correctly) that there is no true 
singularity at these points. The same conclusion holds for the point Zj, for which 
estimated exponent is very close to the positive-integer value A = 2. Finally, it seems 
that z^ = —1 (i.e., q = 0). For this example, at least, the DA method has itself told 
us (correctly) that the only true singularities are those shown in Table 10. 

In Figure 2(a) we plot the limiting curve for m = 3f in the complex 2;-plane 
together with the singularities found by the series analysis. For comparison, we 
have also estimated the radius of convergence of the series from the formula rconv = 
liminf„^oo knl"^''"': despite some oscillations, we are able to obtain the estimate 
Tconv ~ 0.5982, which is close to (and slightly larger than) the modulus {z^l ^ 0.5774 
of the singularity closest to the origin (see Table 10). 



^^In some cases, we find a "broad" peak that has large counts in two (or more) neighboring cells. 
These peaks are counted as one peak; but we have always looked carefully to see whether they have 
a finer structure or not (i.e., if they correspond to two nearby singularities). 
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The case m = 4p is more involved from a numerical point of view because, even 
though the free-energy singularities are of the same type as for m = 3p (namely, 
square- root branch points with A = 1/2), two of them lie extremely close to one 
another. The eigenvalues for this strip are given by the solutions of the quadratic 
equation [6] 

A^-(g^-8g^ + 2V-55g + 46)A + (g^-12g^ + 61g^-169g^ + 269g2-231g + 85) = . 

(5.42) 
The eight endpoints of the corresponding limiting curve are located in the z-plane at 

zi ^ 0.3316354418 ±0.237115247H (5.43a) 

Z2 ^ -0.0681712693 ± 0.4798609413 i (5.43b) 

Z3 ^ 0.2814289723 ±0.4521062477i (5.43c) 

Z4 ^ 0.7398155434 (5.43d) 

Z5 ^ 0.7978491474 (5.43e) 

In Tables 12 and 13 we display our results as a function of A^. (For brevity we 
refrain from showing the histograms.) The three complex singularities {zi, Z2, Z3) and 
their corresponding exponents are well determined with small error bars. However, 
the determination of the two closely-separated real singularities (^4, z^) is (not sur- 
prisingly) rather poor: these singularities can be seen only if we have at least 80 
coefficients, and even with 90 coefficients we obtain only a single singularity located 
in-between the two endpoints. Only for A^ = 100 do we start seeing a (noisy) signal 
of the second endpoint. 

In this case the analysis also gives two singularities that do not correspond to any 
of the endpoints of the limiting curve (see Table 13). One is located aX z = —1 (i.e. 
g = 0) and the other at z ~ 0.285 ± 0.518 i. In both cases, the critical-exponent 
estimate is consistent with A = 0. Therefore, the DA analysis is itself telling us that 
they are not true singularities. 

In Figure 2(b) we plot the limiting curve for m = 4p in the complex z-plane 
together with the singularities found by the series analysis. We have also estimated 
the radius of convergence of the series using rconv = hm inf „_»oo knl"^''"- the result is 
Tconv ~ 0.4216, which is close to (and slightly larger than) the modulus \zi\ ^ 0.4077 
of the singularity closest to the origin (see Table 12). 

5.5.3 Analysis of the series for the bulk, surface and corner free energies 

We are now ready to analyze the large-g series (5.23) ff. for the bulk free energy, 
(5.31) ff. for the surface free energy, and (5.32) ff. for the corner free energy. We 
began by analyzing the series expansions for the exponential of the free energy e-'^*'-^-' 
[i = bulk, surface, corner] using the variable z = l/(g — 1). Later we tried analogous 
analyses for the free energy fi{z) itself, and analyses using 1/q as the variable instead 
of z. 

For the bulk free energy, we have A^ = 47 coefficients. The histograms for the 
K = 1 and K = 2 approximants are displayed in Figure 3(a,b), using the gray-scale 
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coding described earlier. If we compare these histograms to that of the test case of 
Figure 1(a) with a similar number of coefficients (A^ = 50), we see that the zeros do 
not accumulate on small regions on the complex 2;-plane as they did for the test case; 
rather, we find a complex-conjugate pair of broad peaks zi, a broad peak Z2 on the 
real axis, and for i^ = 1 a rather sharp peak 2:3 ^ —1. It is unclear (at least to us) 
what this behavior of the differential approximants is telling us about the singularity 
structure of the function e-'^^""'*^^-'; suffice it to say that we do not expect to obtain 
terribly accurate estimates for the locations and exponents of the singularities. 

For the point Zi using the second-order {K = 2) approximants, our protocol gives 
the estimates 

zi = 0.25(4) ± 0.40(4) i (5.44a) 

Ai = -2.4(6) (5.44b) 

For the first-order {K = 1) approximants, the estimates are 

zi = 0.28(7) ± 0.44(7) i (5.45a) 

Ai = -1.6(6) (5.45b) 

For the broad peak Z2 on the real axis, all the non-defective zeros are precisely 
real (if they weren't, they would come in complex-conjugate pairs so close to the real 
axis that they would fail to be well-separated), and for i^ = 2 our protocol gives the 
estimates 

(5.46a) 
(5.46b) 



For K = 1, the estimates are 



The peak z^, is seen only in the K 



Z2 


= 0.60(8) 


A2 


= -1.5(14) 


Z2 


= 0.64(7) 


A2 


= -1-7(6) 


s K 


= 1 approxir 


^3 


= -1.00(2) 


A3 


= 1-0(1) 



(5.47a) 
(5.47b) 



(5.48a) 
(5.48b) 

The exponent A3 ~ 1 suggests that z^ may not be a singularity of e^^'^^^^^\ but rather 
a simple zero; if so, it would be a logarithmic singularity of the free energy /buik(-2)- 
We also analyzed the series expansion (5.25) of the bulk free energy /buik(-2) [as 
opposed to its exponential (5.23), which we have discussed until now]. In Figure 4 
we show the histogram for the K = 2 approximants of /bulk (-2)- For zi and Z2 the 
picture is roughly similar to what was seen in Figure 3(b); but for z^ we see a much 
sharper peak. Our protocol gives the estimates 

23 = -1.02(2) (5.49a) 

A3 = -0.6(5) (5.49b) 
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The estimate A = —0.6(5) in (5.49b) is barely compatible with the value A = 
corresponding to a logarithmic singularity in /bulk, which is what we would expect if 
g/buik j-^g^s g^ singularity or zero of finite order [e.g. the estimate (5.48b) that suggests 
a simple zero]. 

For the surface free energy, we have N = 47 coefficients (assuming the correctness 
of the conjectural last coefficient). The histograms for K = 1 and K = 2 approxi- 
mants are displayed in Figure 5, and they look even broader than those for the bulk 
free energy. Using second-order [K = 2) differential approximants we are able to 
reliably locate only one singularity, namely 



(The estimate for Z2 is 
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= -1.02(2) 




\2 


= -1.9(15), 



(5.50a) 
(5.50b) 



(5.51a) 
(5.51b) 



which we disregard because the error bar on A2 is so large.) For first-order {K = 1) 
approximants, we find two singularities 

zi = 0.25(5) ± 0.43(9) i (5.52a) 

Z2 = 0.7(1) (5.52b) 

with exponents Ai = —1.3(5) and A2 = —1.6(11), respectively. 

Finally, for the corner free energy, we have A^ = 46 coefficients, and the corre- 
sponding histograms are displayed in Figure 6(a,b). In this case, we find a fairly sharp 
peak around Z2 ~ 0.5 (it is especially sharp for K = 1), and a complex-conjugate 
pair of broad peaks around z ~ 0.2 ± 0.4 i (they are, in fact, slightly less broad than 
for the bulk and surface free energies). As these peaks are close to those found for 
the other two free energies, let us try to locate them carefully. Using second-order 
{K = 2) differential approximants we locate two singularities: 

zi = 0.22(4) ± 0.36(5) i (5.53a) 

Z2 = 0.51(1) (5.53b) 

with associated critical exponents Ai = —1.9(9) and A2 = —1.1(3). These error bars 
are smaller than those found for the bulk and surface free energies (this is especially 
so for Z2), refiecting the sharper peaks in the histogram. For first-order {K = 1) 
approximants, we find the singularities 

zi = 0.23(3) ± 0.36(2) i (5.54a) 

Z2 = 0.509(10) (5.54b) 

with exponents Ai = —1.8(7) and A2 = —1.0(1). As for bulk case, the analysis of 
series for /corner (-2) [rather than its exponential] leads to a new singularity z^ ^ —1, 
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with a quite sharp peak: see Figure 6(c) for the histogram corresponding to the K = 2 
approximants. Our analysis of this peak yields 

^3 = -1.007(8), (5.55) 

with exponent A3 = —0.8(2). The error bar here is about half of that found in 
(5.48a) /(5.49a) for the bulk and surface free energies, reflecting the sharper peak in 
the histogram. 

One possible reason for the better behavior of the corner-free-energy series, at 
least as concerns the real zeros Z2 and Z3, comes from the empirical fact that the 
series for fcomeriq), written in powers of 1/q, has all nonnegative coefficients (as does, 
therefore, the series for its exponential). By contrast, all the other series analyzed 
here have coefficients of both signs, with no clear sign pattern. 

In order to obtain more realistic error bars, we have repeated the analysis using the 
following variations (in addition to the order of the differential approximant K = 1,2): 

• Using the free energy /j instead of its exponential e^\ 

• Using the variable 1/q instead of the variable z = l/(g — 1). 

• "Truncating" the series expansions of e^^^^^ by deleting the terms of order 
z^, . . . , z^ and then dividing by z^^. This truncation ought to have little influ- 
ence on the estimates of the singularities, as the relevant information is encoded 
in the higher-order terms. 

Alas, the estimates obtained with these procedures deviate among themselves by 
about twice the error bars provided by each method. For instance, the position of the 
singularity zi for the bulk free energy takes values ranging from 0.23(3) ± 0.39(3) i 
[for the truncated series expansion of e-^'^""'''^-' with i^ = 2] to 0.30(4) ± 0.51(5) i [for 
the series expansion of fbu\k{^/q) with K = 1]. 

Taking all this information into account, our final estimates for the location of the 
singularities are 

1 = 0.26(4) ± 0.42(4) i 
Bulk ^ ^Z2 = 0.63(10) (5.56a) 

3 = -1.01(2) 

Surface ^ ^ - ^ 0.24(4) ± 0.39(5) . 

^ , = 0.67(14) ^ ^ 

L = 0.23(4) ± 0.37(4) i 
Corner ^ { Z2 = 0.51(1) (5.56c) 

3 = -1.01(1) 

(We refrain from further discussion of the exponents Aj, because the error bars are so 
large.) 
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We see from (5.56) that the estimated locations of the singularities associated to 
the three free energies are (as expected) compatible within errors. Combining these 
estimates, we obtain 

zi = 0.24(4) ± 0.39(4) i [q^ = 2.14(19) ± 1.86(19) i] (5.57a) 

Z2 = 0.51(2) [ga = 2.96(8)] (5.57b) 

Z3 = -1.01(1) [gs = 0.01(1)] (5.57c) 

The absolute values are therefore \zi\ = 0.46(5), \z2\ = 0.51(2) and {z^l = 1.01(1). 

In Figure 7(a,b,c) we compare the estimates (5.56) for zi, Z2, z^ with the two best 
currently available approximations to the square-lattice limiting curve Boo (namely, 
the curves Bm coming from m = 7 with toroidal boundary conditions [10] and m = 11 
with cylindrical boundary conditions [7]). We also show a circle with radius Tconv ~ 
0.4290 [resp. rconv ~ 0.4225, rconv ~ 0.4139] indicating the radius of convergence 
obtained by estimating rconv = liniinf„^oo |on|~^^"- 

Let us discuss the three singularities in order. 

The singularity of smallest modulus is the complex point zi, given by (5.57a). 
From Figure 7 we see that the estimated location zi lies fairly close to the most 
prominent T point of the limiting curve Bi^^^ (solid black line), namely [10, sec. 3.6] 

^T,7,„^ = 0.26158090412±0.3162213774i [gT,7to. = 2.5531414480±1.8775702667i] , 

(5.58) 
and one can conjecture that this is not an accident. On the other hand, this singularity 
lies rather far from the T point of the limiting curve Bll^^^ (dashed green curve), 
namely [7, sec. 3.3] 

ZT,ii,,, = 0.409 ±0.218 2 [gT,ii,,, = 2.902 ± 1.015 2] . (5.59) 

Thus, we are unable to draw any firm conclusion about the physical nature of this 
singularity. It does, however, seem likely that this singularity lies somewhere on the 
limiting curve B^o — possibly (but not necessarily) at a "special" point of B^o such 
as a T point. 

The two singularities on the real axis are given by (5.57b,c). The estimate for Z2 
is compatible within errors with q2 = 3, which is the physical critical point gc(sq) = 3 
of the square-lattice zero-temperature Potts antiferromagnet [68-75] and is expected 
to be the uppermost real crossing point of the limiting curve B^o (see Section 6.2 
below). The estimate for z^ is compatible with ^3 = 0, which is expected to be the 
lowermost real crossing point of the limiting curve Boo- 

In conclusion, the analysis of the free-energy series expansions has led to the loca- 
tion of four singularities in the complex g-plane: two of them are the real points q = 
and q = 3, which correspond to the transitions between the disordered and critical 
phases, while the other two are a pair of complex-conjugate points that presumably 
lie somewhere on the limiting curve Boo, possibly at a "special" point such as a T 
point. 
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Our difficulties in locating the singularities of the free energies — in particular, 
the wide dispersion of estimates visible in the histograms (Figures 3-6), which con- 
trasts strongly with the behavior observed in our test cases 3p and 4p (Figure 1) 
— suggest that the singularity structure of the infinite-volume free energies may be 
quite different from that of the finite-width strips that we used as test cases. In 
particular, it is possible that the entire limiting curve Boo is a natural boundary for 
the infinite-volume free energies — possibly with a very soft singularity — and that 
the differential approximants are trying, with limited success, to mimic this behavior. 
In this interpretation, the points Zi, Z2 and z^ would simply be approximations to 
the places on Boo where the free energies are most strongly singular; and one would 
furthermore expect that as the number A^ of series coefficients increases, the number 
of peaks in the histograms will increase, and the zeros will gradually condense on the 
whole of the curve Boo (perhaps with high density on some parts of the curve and 
low density on others). It would be interesting to try further numerical tests of this 
scenario; but they will probably require fairly radical extensions of our large-g series, 
e.g. to order N ^ 100. 

6 Numerical results for widths m = Of, IOf, Hf? 12f 

As part of this work, we have computed the transfer matrices for square-lattice 
strips with free boundary conditions and widths m = 9, 10, 11, 12. (Widths m < 8 
were computed in a previous paper [6].) In this section we use these newly-computed 
transfer matrices to study the real and complex roots of the chromatic polynomials 
Pmpxnp{(l) for 771 = 9, 10, 11, 12, focusiug on the behavior in the infinite-length limit 
(n — s> cxd). In Section 6.1 we discuss the limiting curve Bm and the isolated limiting 
points for each width m. In Section 6.2 we study the behavior of the real crossing 
points qoi'm) and attempt to extract q'o(sq) = lim qo{m). In Section 6.3 we discuss 

the real chromatic roots for finite lattices m-p x n-p. 

6.1 Limiting curves and isolated limiting points 

Having computed the transfer matrices for m = 9f, lOp, Up, 12f, we can analyze 
them as in [6-10] to extract the chromatic roots for finite-length lattices (n < oo) 
as well as the limiting curves and isolated limiting points in the infinite-length limit 
(n -^ cxo). Unfortunately, we have not been able to compute the full limiting curves 
Bm, as this would have required a major computational effort. However, we have tried 
to locate the most important points on each of them, i.e. crossings of the real axis 
and endpoints near the real axis, using the direct-search method [6]. These results 
are summarized in Table 14.'^'^ (Furthermore, the chromatic roots for n = bm and 
n = 10m shown in Figures 8(a)-(d) give a fairly good idea of the general shape of 



^Opor widths m < 11, wc implemented the direct-search method in Mathematica, using 100-digit 
internal precision. For width to = 12, we used the double-precision Fortran subroutines in the 
ARPACK package [76]. We also checked our results for 7 < to, < 11 using ARPACK; the results are 
consistent with those obtained using Mathematica but are less precise. 
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the curve Bm-) Finally, for all these strips there are real isolated limiting points at 
the Beraha numbers q = 0,1,2, B^, and for m = 9,11,12 there also appear to be 
complex isolated limiting points (we cannot, however, guarantee that we have found 
all of them). In Figures 8(a)-(d) these points are marked with a x sign. 

For m = 9, the dimension of the transfer matrix is 179. In Figure 8(a) we have 
plotted the chromatic zeros in the complex g-plane for the strips 9f x 45f and 9f x 90f- 
The limiting curve crosses the real g-axis at a single point: 

go(9) ^ 2.70165995678. (6.1) 

Figure 8(a) also suggests that there is a complex-conjugate pair of isolated limiting 
points at g ^ 2.5946 ± 0.5963 i. 

For m = 10, the dimension of the transfer matrix is 435. In Figure 8(b) we have 
plotted the chromatic zeros for the strips lOp x 50f and lOp x IOOf- In this case 
the limiting curve does not cross the real g-axis; rather, there is a pair of complex- 
conjugate endpoints very close to the real g-axis, at 

go(lO) ^ 2.7343903604 ±0.0003924978 2. (6.2) 

For m = 11, the dimension of the transfer matrix is 1142. In Figure 8(c) we have 
plotted the chromatic zeros for the strips lip x 55f and lip x llOp- The limiting 
curve crosses the real g-axis at 

go(ll) ^ 2.7608973951 . (6.3) 

Figure 8(c) also suggests that there are two complex-conjugate pairs of isolated lim- 
iting points at g ^ 2.6555 ± 0.4978 i and g ^ 2.4648 ± 1.1380?. 

Finally, for m = 12, the dimension of the transfer matrix is 2947. In Figure 8(d) 
we have plotted the chromatic zeros for the strips 12f x 60f and 12f x 120f. In this 
case the limiting curve does not seem to cross the real g-axis; rather, there is a pair 
of complex-conjugate endpoints very close to the real g-axis, at 

go(12) ^ 2.782817590 ± 0.00018700 i . (6.4) 

Figure 8(d) also suggests that there is a complex-conjugate pair of isolated limiting 
points at g ?« 2.6172 ± 0.7562 i 

In view of the results reported in Table 14, we conjecture that the limiting curve 
Brn crosses the real axis for all odd m and has a complex-conjugate pair of endpoints 
very near the real axis for all even m > 8. 

6.2 Value of go(sq) 

We can try to use the results reported in Table 14 to obtain the value of go(sq) = 
lim go(m). Of course, we expect go(sq) = 3 [68-75], but it will be interesting to see 

with what accuracy this result can be obtained. 



Let us also warn the reader that the value for to = 8 reported in [6] is wrong; the correct value is 
displayed here in Table 14. 
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As just discussed, there are clear parity effects in go(^): the hmiting curve Bm 
crosses the real g-axis for odd m but not for even m. We have therefore split the data 
into two sets according to the parity of m, and analyzed each set separately. Please 
note that the data for qoiiTi) are essentially exact: there is a tiny non- statistical error 
of order 10~^° in their numerical estimates. In order to keep to the standard notation 
in finite-size-scaling theory, we here denote the strip width by L instead of m. 

For each data set we have considered several Ansatze, and for a given Ansatz 
with k free parameters, we have taken into account the points with L = Lmm, -^min + 
2, . . . , Linin + 2(/i; — 1). So in each fit there are no degrees of freedom. From the 
variation of the estimates as Lmin is increase, we can roughly estimate the error bar. 

If we use the power-law Ansatz 

qo{L) = qo{sq)+BL-^ (6.5) 

for the data coming from odd L, we obtain the following estimates using Lmin = 7 

go(sq) = 2.999(5) (6.6a) 

A = 1.108(16) (6.6b) 

where the error bars are defined to be twice the distance between these estimates and 
those obtained with Lmin = 5. 

If we play the same game for the quantity Reg for even L, we arrive at the 
estimates for Lmin = 8 

go(sq) = 3.00(7) (6.7a) 

A = 1.10(29) (6.7b) 

where again the error bars are twice the difference from the estimates for Lmin = 6. 
The convergence in this case is not very good, so we cannot trust these results too 
much. 

The results (6.6)/(6.7) are consistent with the behavior 

go(sq) = 3 (6.8a) 

A = 1 (6.8b) 

Of course, we expect that go(sq) = gc(sq) = 3; moreover, the correction exponent 
A = 1 is to be expected for free b.c. because of the surface effect. For cylindrical 
boundary conditions [7, Section 4], the estimates are less well converged but are 
consistent with go(sq) = 3 and suggest that A > 1.7 (e.g. perhaps A = 2). Finally, 
very recent work [10, Section 5.1] gives strong evidence that go(sq) = 3 holds also for 
toroidal boundary conditions. ^^ 



^-'^The authors of Ref. [10] found that for toroidal square-lattice strips of odd widths 3 < L < 11, 
the value oiqo(L) is exactly equal to the expected limiting value (7o(sq) = 3, while for toroidal strips of 
even widths 2 < L < 12, the fits of the numerical data to the Ansatz (6.5) gave (/olsq) — 2.999±0.012 
and A = 2.04 ± 0.08, in agreement with the conjectured behavior (6.8a) but suggesting that here 
A = 2. 

59 



Better fits can in principle be obtained by imposing the values (6.8). In particular, 
we expect the following Ansatz to describe better the data points: 

go(L) = 3-BL'^-CL~^. (6.9) 

However the results are not very stable: for odd L we find that for Lmin = 7, -B ~ 
2.436, and C ^ 2.243, while for Lmin = 9, B ^ 2.283, and C ^ 2.719. A similar 
behavior is found for even L: for Lmin = 8, -B ~ 2.409, and C ~ 2.473, while for 
Lmin = 10, -B ~ 2.357, and C ~ 2.994. We conclude that our numerical data are 
not accurate enough to be able to accurately determine corrections to scaling to the 
behavior of qo{L). 

For odd L, we can also look at Imgo(-^)- In R-ef. [7, Conjecture 4.2] it was conjec- 
tured that for square-lattice strips with cylindrical boundary conditions Imgo(-^) ~ 
i?5 = r^^ where r = (1 + v^)/2 is the golden ratio. We have fitted our numer- 
ical results for Imgo(-^) with L = 8,10,12 to check whether this conjecture holds 
or not also for free boundary conditions. A 3-parameter fit to a power-law Ansatz 
Imgo(-^) = A + BL~^ shows that A ^ —0.00042 is indeed a small number. We 
also obtain that B ^ 0.033 and A ^ 1.606. As we have only three data points, 
there are no degrees of freedom in the fit and the error bars cannot be reliably es- 
timated. In conclusion, we expect that Imgo(-^) decays to zero exponentially fast. 
Thus, we use the improved Ansatz: Imgo(-^) = Ax B^. The results for Lmin = 10 are 
A = 0.16(13) and B = 0.690(75). The error bars are computed as twice the difference 
between the estimates for Lmin = 10 and Lmin = 8. If we compare B~^ = 1.44(16) to 
r ^ 1.6180339887, we see that they are barely compatible within errors. As the error 
bar for B~^ is 20 times larger than the corresponding error for cylindrical boundary 
conditions [7, Section 4], we cannot draw any firm conclusion about the validity of 
the above conjecture for square-lattice strips with free boundary conditions. 

6.3 Real chromatic roots 

It is also of interest to study the real chromatic roots for lattices of finite length 
n. For every length n there are, of course, roots at g = and q = 1] and there are 
also roots converging (exponentially rapidly) as n — >■ oo to the isolated limiting point 
at g = 2. Here we shall concentrate on the roots converging to the isolated limiting 
point at g = -B5 = (3 + \/5)/2 and, for m = 9 and 11, to the real crossing point qo{m). 

For width m = 10, 12, the real roots converging to Br, do so monotonically from 
below; we refrain from presenting the details. There are no real roots above B^. 

For widths m = 9 and 11, by contrast, two interesting things happen (see Tables 15 
and 16). On the one hand, the real roots converging to -B5 do so with parity (— l)""*"^, 
i.e. alternating from above and below (starting from n = 39 and n = 23, respectively). 
On the other hand, for odd lengths n, there are also real roots converging to qo{iTi) 
from below. We conjecture that these behaviors will persist for all larger odd widths 
m. 

These results for m = 9, 11 provide a (presumably infinite) family of counterex- 
amples to a conjecture made in [6] concerning the real chromatic roots of bipartite 
planar graphs. For further discussion and for a revised conjecture, see Appendix B. 
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An important question is the rate of convergence of these real zeros towards the 
hmiting values -B5 and qo{m). For B^ we expect an exponentially rapid convergence 
(as with all isolated limiting points in the Beraha-Kahane- Weiss theorem); for qQ[m), 
by contrast, we expect a 1/n convergence. We have tested these predictions by fitting 
the given real roots g„ to a power-law Ansatz 

q = q^ + An-^ , (6.10) 

using the data points n > a variable threshold n^nm- For the sequences of zeros con- 
verging to B5, the estimates of the power A appear to increase without bound as nmin 
is increased (e.g. A > 45), suggesting that the convergence is indeed exponentially 
fast. For the sequences converging to qo{m), by contrast, we obtain powers A ^ 1.1 
for small values of n^^i^, which slowly decrease toward 1 as we increase nmin. In par- 
ticular, for nmin = 99 we obtain A ^ 1.023 for m = 9 and A ^ 1.037 for m = 11. 
This is consistent with the 1/n behavior expected from the Beraha-Kahane- Weiss 
theorem. 



7 Summary and open problems 

In this paper we have derived some new structural results for the transfer matrix 
of square-lattice Potts models with free and cylindrical boundary conditions. In 
particular, we have obtained explicit closed-form expressions for the dominant (at 
large |g|) diagonal entry in the transfer matrix, for arbitrary widths m, as the solution 
of a special one-dimensional polymer model. We have also obtained the first 47 (resp. 
46) terms in the large-g expansion of the bulk and surface (resp. corner) free energies 
for the zero-temperature antiferromagnet (= chromatic polynomial). Finally, we have 
computed chromatic roots for strips of width m = 9, 10, 11, 12 with free boundary 
conditions and located roughly the limiting curves. 

Let us end by listing some open problems for future research. The first group 
of such problems would be to understand theoretically (and ideally to prove rigor- 
ously) some of the regularities that were observed empirically in Sections 4 and 5, for 
instance: 

1) The empirically observed property that the polynomials p^{k) and p^{k) arising 
in (3.42)/(3.71) ff. have for even s > 2 a factor k — s/2. 

2) The formulae (3.43b) and (3.72a,b) relating the coefficients a^{k + s) and d^{k + 
s) for s = 1 and s = 0, 1, respectively, to Fibonacci numbers (and possible extensions 
of these formulae to higher values of s). 

3) The fact that a^^ and a^^ [defined in (3.48) and (3.75)] are, for each fixed i > 0, 
the restriction to integers A; > £ of a polynomial in k of degree £. 

P /P 

4) The behavior (4.6)/(4.7) relating the leading eigenvalue A* (m) to the domi- 
nant diagonal entry tp/p{m). 

5) The fact that the coefficients cj (m) in the large-g expansion of logA^ (m) 

P /P 

[cf. (4.9)] are the restriction to integers m > m^^^^ [k) of a polynomial of degree 1 in m. 
Finite-size-scaling theory (Section 5.1) gives a nonrigorous (but physically intuitive) 
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explanation of why the c^ (m) are polynomials of degree 1 in m for large enough m. 
But we do not really understand why the cutoffs rnj^:^^ {k) take the values they do; 
even less do we understand the improvements found empirically in Sections 4.2 and 
4.3. And of course we lack a rigorous proof for all of this. 

6) The correction term of order ^^Lmax-i found empirically for the surface free 
energy [cf. the paragraph preceding (5.31)]. 

Among the more fundamental extensions of this work, we would like to propose 
the following: 

1. Obtain longer large-q series. It is natural to ask whether our large-g series 
for the bulk, surface and corner free energies (Section 5.2) can be extended to higher 
order. Modest extensions are no doubt possible over the next few years, with in- 
creased computer power; but the memory and CPU-time requirements are heavy, 
even for width m = 12 or 13. It seems to us that radical algorithmic improvements 
in the implementation of the transfer-matrix method will be required in order to go 
to significantly larger widths m. Unfortunately, our results from the series analysis 
(Section 5.5.3) suggest that these series are badly behaved, and that it will be nec- 
essary to obtain many more terms — for instance, doubling the length of our series 
from N = 4:7 to N K, 100 — in order to improve significantly our understanding of 
the analytic structure. 

2. Extend this work to the triangular lattice. Extension of this work to the trian- 
gular lattice is particularly interesting in view of Baxter's [77, 78] conjectured exact 
solution for the bulk free energy of the triangular-lattice chromatic polynomial (see 
also [8, Section 6] for a critical discussion). It ought to be fairly easy to obtain exact 
formulae for the dominant transfer-matrix entry, analogous to those obtained here 
in Section 3 for the square lattice. Likewise, the transfer-matrix calculations in Sec- 
tions 4 and 5 can easily be extended to the triangular lattice, as discussed in [8], 
albeit possibly for slightly smaller widths. The main trouble arises in the use of the 
finite-lattice method: using transfer matrices of widths up to Lmaxj we obtain the 
large-g series through order 4Lmax — 2 on the square lattice (Section 5.2) but only 
order ^ 2Lmax on the triangular lattice [67, Section 3.1]. 

But the calculation of the large-g series for the hulk free energy is expected to be 
superfluous, as Baxter [77, 78] has an exact expression 

x-Ll (i_a;6j-5)(i_a.6i-4)(i_2,6j)(i_a.6j+i) y'- ) 

(where g = 2 — x — x~^ and |a;| < 1) that is claimed to represent the exponential 
of the bulk free energy at large enough |g| (namely, outside the limiting curve Boo). 
The finite-lattice calculations would thus serve only to check (7.1) through order x^'^^ 
(and also to compute the series for the surface and corner free energies). Assuming 
that this check confirms (7.1) — as we expect that it wilF^ — we will learn nothing 



^^See also the numerical confirmations in [8, Section 6.3]. 
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new! One could equally well carry out a different ial-approximant analysis like that of 
Section 5.5.3 directly on the series obtained by expanding (7.1), for which it is easy 
to obtain 1000 or more terms. 

On the other hand, there is little point in carrying out such a differential-approximant 
analysis, since we already know the analytic structure of (7.1): namely, it is analytic 
and nonvanishing in the disc |a;| < 1, or equivalently in the complex g-plane with a 
cut along the real interval < g < 4. It follows that, if Baxter's formula (7.1) is 
correct at large l^l, then the limiting curve B^o for the triangular-lattice chromatic 
polynomial is not a natural boundary for the bulk free energy — contrary to what we 
suspect to be the case for the square-lattice chromatic polynomial (see the discussion 
at the end of Section 5.5.3). 

When all is said and done, our lack of understanding of the triangular-lattice 
chromatic polynomial concerns small \q\ — namely, the region inside the limiting 
curve Boo — and not large |g| (see [8, Section 6] for discussion). This issue is probably 
best addressed by trying to understand better the structure of Baxter's [77, 78] Bethe- 
Ansatz solution for finite widths m, including the effect of boundary conditions. 

3. Carry out transfer-matrix calculations for the diced lattice. The discussion 
in Appendix B and in [79] suggests that it would be very interesting to carry out 
transfer-matrix calculations analogous to those of [6-10] and the present paper for 
the diced lattice (namely, the dual of the kagome lattice). We have established in [79] 
that gc(diced) > 3, but we do not know whether the limiting curves Bm intersect the 
real axis at some value qoitn) > 3, or merely have complex-conjugate endpoints that 
tend to gc(diced) as m -^ oo. 

4. Understand the analytic nature of complex phase boundaries. As discussed at 
the beginning of Section 5.5, it is proven in some cases, and expected in general for 
statistical- mechanical models in dimension d > 2, that the infinite- volume free energy 
has a "soft" essential singularity (where it is infinitely differentiable but not analytic) 
whenever a first-order-phase-transition point in the physical region is approached. 
But it is unclear whether the complex phase boundaries are natural boundaries of 
the free energies defined in each phase. (In our case, we are uncertain whether the 
limiting curve B^o is a natural boundary for the bulk free energy of the square-lattice 
chromatic polynomial.) It would be useful to resolve this question in at least one 
model (for instance, the analyticity or nonanalyticity at pure-imaginary magnetic 
field of the low-temperature Ising ferromagnet), ideally by a mathematically rigorous 
proof. 

A Some combinatorial identities 

In this appendix we will prove Propositions 3.5 and 3.7 concerning the quantities 
al{k — s) and a^{k — s), respectively. We will also prove eq. (3.43a) for the quantity 
alik). 

We shall use Knuth's notation for falling powers [42] 

x^ = x{x-l)---{x-n + l) . (A.l) 
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We also use the standard convention [42] for the definition of binomial coefficients: 

(A.2) 



X \ I x-/k\ for integer A; > 



for integer /c < 



where x can be a real or complex number (or more generally an algebraic inde- 
terminate) while k is always an integer. Finally, we adopt the convention that 
k\ = T{k + 1) = oo for integer k < 0, so that (A.2) could be written simply as 

Q=x^/ki 

A.l Some preliminary lemmas 

Let us start with a simple lemma that will make the analysis easier: 

Lemma A.l For any integers s,a>0, let us define the polynomial 
Then 



r\ [s — r)\ 

r=0 ^ ' 



(hence in particular Ts^a{^) = for a > s). 

Proof. Expanding r- in (A. 3) and using the fact that r- = for a > r, we obtain 

■^ — ' r\ [s — r]\ 

r=a ^ ' 

If a > s, then the sum is empty, so let us suppose that a < s. Then we can rearrange 
the factorials to obtain 



z. 


a 


(^ 


r-a)\{s- 


-r)! 


(^ 


1 


a)\ 


FA'- 

r=a ^ 


-fly 


(^ 


A" 
A" 


a)\ 


r'=0 ^ 

(A + 1)^-" 


— a 
r' 



(A.6) 



D 



The next lemma will allow us to perform the inner sums in (3.34) and (3.63): 
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Lemma A. 2 For any integers i,a,b > 0, let us define the polynomial 



r=0 

Then 



^...^±frr"-%-r) 



(A7) 



^'--W = £!(!"+ i-QI ^'-'"-^'- ^^' 



where 

min{a/) 



Q^^(6,A) = y ( "" ) (a + 6-£)^i^£^A™(A-l)^-™ (A.9) 

is a 'polynomial in b, A that is of degree a in the variable b and of degree i in \, with 
leading coefficient 

mQe,a{bA) = (A -if (A.IO) 

and hence [b"'X^]Qe,a{b,X) = 1- 

Proof. Let us first consider the case i — a — b — 1 > 0. Then ( ~''~ ~ ) = whenever 
r > i — a — b—l. On the other hand, (^_ J = whenever i — r>b, i.e. r < £ — b. But 
these two conditions on r involve all possible cases, as a > 0. Therefore, Fia^bi^) = 
for £ - a - 6 - 1 > 0. 

Let us now consider the nontrivial case i — a — b— 1<0. As the upper index in 
( ~"r ~ ) ^^ ^'^^ negative, it is convenient to rewrite (A. 7) as 

v^ , . •, ^ ff^ — ^ + b + a\ / b \ 



r=0 



where now all the indexes are nonnegative. We can then rewrite the binomial coeffi- 
cients in terms of factorials as follows: 

^ Y^. y {r-i + b + a)\ b\ 

7^0 r\{a + b-t)\ (£-r)!(6-£ + r)! 

The product in (A. 12) can be written as a sum using the "binomial theorem" for 
falling powers [42, Exercise 5.37]. Defining a; = 6 — £, we have 

a 

JJ(r + a; + m) = [r + x + a)- (A.13a) 

^(''\r^{x + a)'^=^ (A. 13b) 



m=l 



m=0 
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and therefore 



^' ^r«^.„,. Ma-™v-Mr^^ 



Q,,,(6,A), (A.14c) 



m=0 

6! (-1 



{a + b-i)\ i\ 



where we have used Lemma A.l and the definition (A. 9) of Qe,aib,X) [note that 
the sum in (A. 9) can be stopped at m = £ because i— = ior m > i]. Clearly, 
Qe,a{b, A) is a polynomial in b, A that is of degree at most a in the variable b and of 
degree at most i in A. Indeed, this is its exact degree, as the leading coefficient is 
[b'']Qi,a{b, A) = (A - 1)^ [from the m = term in (A.9)] and hence [b^\^]Qi,a{b, A) = 1. 
D 



Remarks. 1. When A = 1, (A. 7) reduces to the Vandermonde convolution 



while (A.9) reduces to 



Fi,a,i.a) = r ! M , (A.15) 



a\ . . r ..a-e .. «' (O' + b 



QeAb, I) = [,]ia + b-ir-^i\ = ^--^ ^ '^, " , (A.16) 



so that (A. 8) gives 



Fi,a,a) = i-iy {''^ = (^ ! ^) ■ (A.17) 



2. In our applications we will always take i = k — p > and A = 3. In addition, 
in Theorem A. 4 we will take (a, b) = (s — l,p + 1) with p > 0, and (a, b) = {s,p — 1) 
with p > 1; finally in Theorem A. 6, we will take (a, 6) = (s — l,p — 1) with p > 1. 

A. 2 Proof of formula for a^{k — s) 

The next lemma deals with the outer sums we find in (3.34): 

Lemma A. 3 Let a = or 1. Then, for any integers s > 1 — a and k > s + a, let us 
define 

( — 1)^ ] Fk-p,s+a-l,p+l-2a{^) , (A. 18) 



p=a 
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where Fi^a,b{^) is defined by (A. 7). We then have 



E ^ ^lu^.^L '^ Q.s^iik+l-q,X) ^fa = 



A.a,s,kW — ■* 



,=0 1' (" - «)' 



(A.19) 



if a = 1 



where the polynomial Qe,a{b,X) is defined by (A. 9). In particular, Aa,s,fc(A) is the 
restriction to integers k > s + a of a polynomial in k, A that is of degree s in k and 
A separately, with leading coefficient 



[A;1A«,.,fc(A) 



-1)^ 



S! 



(2-Ar 



(A.20) 



and hence [k"^ \^] Aa,s,k{X) = (— l)*^"/s!. The next subleading coefficients are given 
by 



[k'-']Aa,s,kW 



(2 — A)''^^ I (4 — A — 2s + As) if a 







2fs-l)! 



(2-A)(s + l) ifa = l 



(A.21a) 



[k-'-'] A„,.,fc(A) 



(2 _ \y-i [ -(3As2 - 6s2 - 7As + 2s + 2A - 4) if a = 



24rs-2)! 



[k'-'] A„,s,fc(A) 



(2 - A)^-^ 
48 (s- 3)! 



X 



3As2 - 6s2 + 5As - 34s + 2A - 4 if a = I 

(A.21b) 

-48 + 80s - lOAs - 3A2s - 20s2 
+2As2 + 4A2s2 - 4s3 + 4As3 - A^s^) if a = 



[A;^-1 A,,s,fc(A) 



(2 -A) 



(2 - A)s(s + 1)(26 - A + 2s - As) 



' (2 - A) (-16 + 8A + 924s + 18As 
-1310s2-125As2 + 420s3 
+90As3 + 30s4 - ISAs^) 



s-2 



if a = 1 
(A.21c) 



if a = 



5760 (s- 4)! 



X < 



32 - 32A + 8A2 + 1032s + 888As 
+18A2s - 5300s2 + 1220As2 
-SAls^ _ I560s^ + 840As^ 
-SOX'^s^ - 60s^ + 60As^ - ISA^s^ if a = 1 

(A.21d) 
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Proof. In (A. 18) let us insert the result (A. 8) of Lemma A. 2 for Fk-p,s+a-i,p+i--2a{^) ',^^ 
we get 



^a,s,fc(A) 



k 
p=a 



- S 
P 

xQ 



p — 1 
a 



-lf{p+l-2a) 



{2p + s 

k-p,s+a-l{P + 1 



k-a)\{k-p)\ 
2a, X) . 



(A.22) 



First of all, we notice that all terms with p < k — s — 1 vanish: for ifO < k — s—p — 1 < 
p-a, then (''"^1^"^) = 0, while if 2p + s - A; - a < 0, then l/{2p + s - k - a)\ = 0; 
and these two cases cover all values oi p satisfying p < k — s — 1. So the only nonzero 
contributions to the sum (A.22) come from p > k — s > a. Using the variable 
q = k — p, we can rewrite (A.22) as 



Aa,s,k{X) 



E 

q=0 



g=0 



g=0 



q-s-l\ (-1)'= {k + l-q- 2a)! 



k — q — aj {k + s — 2q — a)\ q\ 
X Qq,s+a-i{k + 1 - g - 2a, A) 
k + s-2q-a\ (-1)'?+" {k + 1 - q - 2a)! 



k — q~a J {k + s — 2q — a)\q\ 
X Qg^s+a-i{k + 1 - g - 2a, A) 

-1)"+" (A;-g-2a + l)! 



q\ {s — q)\ (k — q — a)\ 



X Qq,s+a-i{k + 1 - g - 2a, A) 



(A.23) 



Evaluating this for a = 0, 1 yields (A. 18). This expression for Aa,s,k{X) is manifestly a 
polynomial in k and A of degree at most s in each of the variables separately. Indeed, 
using (A. 10) and performing a binomial sum over q, we obtain (A. 20). 

The computation of the subleading terms (A. 21) is more involved. Let us start 
with the case a = 0. Our goal is to extract the leading powers of k in (A. 19), so 
our first task is to convert the falling power {k + s — 2q) '^~^~"^ that appears in the 
definition of Qg^s-i{k + 1 — g. A) [cf. (A. 9)] to a regular power using the Stirhng cycle 
numbers [42]. The key formula is obtained from (3.47) and the binomial theorem: 



(k + x) 



p=0 



r=0 



m 
m — p 



m-p\ 

X ^ . 
r — p ' 



(A.24) 



■^■^Please note that the hypotheses of Lemma A. 3 {a S {0, 1}, s > 1 — a and k > s + a) together 
with the summation hmits a < p < k imply that the three subscripts in i^fe-p,s+Q-i,p+i-2a are 
nonncgative; this justifies the use of Lemma A. 2. 
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Using this formula in the definition of Qq^s-i{k + 1 — g, A) allows us to obtain the first 
five leading terms in the expansion of the summand in (A. 19): 

{k + l-q)Qg^s-iik + l-q,\) 

= (A; + 1 - g) ^ ("^ " ^Va; + s - 2g)^=i^g^A'" (A - 1^"^ 

4 

= J]F-'af) + 0(F-^) , (A.25) 
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where the coefficients a\ are given by 

4°) = {\-iy 



M 



,(0) 



(A-1) 



9-1 



-q (1 - 2s + As) + -L-J 5(5 + 1) 



(A.26a) 
(A.26b) 



(A-1) 



q-2 



s-1 



g^-^^-^(2-2s + As) 



-Q 



(A -2) 
2 

(A-1)^ 
24 



(-2 + 3s - 2s=^ + As^) 
s(s + l)(2 + 3s) 



(A.26c) 



,(0) 



(A-1) 



+q- 



9-3 



s-1 



_g3ii_J)!(3_25 + As) 



(A-i: 



16 - 6A - 16s + lOAs - A^s + 4s^ - 4As^ + A^s 



2 I \2„2\ 



-q ^ ^ (12 - 25s + 17s2 - As^ - 6s=' + SAs^) 

(A-l)\2,. , „^2 



24 



s^ 1 + s 



(A.26d) 



,(0) 



(A-1) 



9-4 



s-1 



giIi^(4-2s + As) 



_^3 (A l)jA 2) ^^^ _ ^2^ _ 2g^ ^ ^^^^ _ ^^2^ ^ 4^2 _ ^^^2 ^ ;^2^2^ 



+r 



2 

(A -1)2 



-208 + 56As + 236s - 102As + 4A2s - SSs^ + SSAs^ 



-TA^s^ + 12s^ - 12As^ + 3A2s^) 
^^ ^ (s - 1)(8 - 14s + 7s2 - 2s3 + As=^ 



-q 



+ ^-7t::7^s(s + l)(-8 - 10s + ISs^ + ISs^) 



960 



(A.26e) 
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In the derivation of these formulae we have used the following special values of the 
Stirling cycle numbers [80, Chapter 5, Exercise 16] (see also [30,42]): 



n 
n — 1 

n 
n-2 

n 
n — 3 

n 
n-4 



n\ 




2J 




n^ 


3n-l 


sy 


4 


n^ 


fn\ 


4y 


V2; 


n^ 


ISn^ - SOn^ + 5n + 2 


5) 


48 



(A.27a) 
(A.27b) 
(A.27c) 
(A.27d) 



The second task is to perform the sum over q in (A. 19). It is clear from (A. 26) that 
the coefficients a\ are polynomials in q; and we have chosen to write them in terms 
of the falling powers q— in order to facilitate the use of Lemma A.l to perform the 
sum over q in (A. 19). Each coefficient al in (A. 26) contributes only to the coefficient 
[k'-']Ao,s,kW: 



[k'-']Ao,s,kW = J2 



-lY 



q=0 



g!(s-g)! 



a,- 



(0) 



for i > 0. 



Using the definition (A. 3) of Tsa(A), we can rewrite (A. 26) as 
[F]Ao,,,fe(A) = T,,o(l-A) 



(A.28) 



(A.29a) 



[k'~']Ao,s,k{X) = (A-1)-^ 



;i-2s + As)T,,i(l-A) 



+^^ir^s(s + l)T,,o(l-A) 



(A.29b) 

along with similar formulae for [A;''~*]Ao,s,a:(A) with i = 2,3,4. From these equations 
together with Lemma A.l we can derive after some algebra the results (A.20)/(A.21) 
for a = 0. 

The computation for a = 1 is similar. The five leading terms for Qq^s{k — 1 — g. A) 
can be written as 

Q,,,(A;-l-g,A) = ^ T^Va; + s - 2g - l)^^g^A™ (A - 1)^"™ 

m=0 ^ ^ 

4 

= ^F-*af^ + 0{k'-^) , (A.30) 



i=0 
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where the coefficients a\ are given by 



M) 



.W 



(A-1) 



9-1 



(A-1) 



g-2 



q(2-\)+^-^—^{s-l) 



q^ (2 -\f -q{\- 1)(6 - 2A - 2s + As) 



(A.31a) 
(A.31b) 



+ ^^TTT^(S-2)(3S-1) 



12 



(A.31c) 



,(1) 



(A - 1)5-3 



s 
3 
(A-1) 



q^{2-Xy -q 



3(2-A)(A-l) 



10-3A-2s + As) 



4 
(A -1)3 



[s-3)(20-4A-6s + 3As) 



s(s-l)(s-3) 



(A.31d) 



a? 



(A-1) 



g-4 



g^ (2 - A)^ - 2q^ (2 - A)2(A - 1) (14 - 4A - 2s + As) 



+q^- — r-^ (328 - 208A + 28A2 - 124s + lOOAs - lOA^s 

+12s2 - 12As2 + 3A^s2) 
(A-p3 






2 

(A-1)^ 



240 



(s - 4)(-18 + 2A + 12s - 3As - 2s2 + As^) 
s-4)(2 + 5s-30s2 + 15s3) 



(A.31e) 



We now perform the sum over q in (A. 19). Again, each coefficient a\ in (A. 31' 
contributes only to the coefficient [k^''^]Ai^s,k{^)'- 



D" 



[knAi,s,kW = -Y. ^u.-nv ^^'' f^^^^o. 



q=0 



g!(s-g)! 



(A.32) 



Using the definition (A. 3) of T5^a(A), we can rewrite (A. 31) as 
[F]Ai,,,fc(A) = -T,,o(l-A) 



(A.33a) 



[F-VM,fe(A) = -(A-l)'^ 



(2-A)T,,i(l-A) 



+^^i7^(s-l)T.,o(l-A) 



(A.33b) 
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along with similar formulae for [k^ '^]Ai^s,k{^) with i = 2,3,4. These equations lead 
via Lemma A.l to (A.20)/(A.21) for a = 1. D 

We now proceed to prove the main result of this subsection: 
Theorem A. 4 (= Proposition 3.5) For each integer s > 1, the quantity 

is, when restricted to k > s + 1, given by a polynomial in k of degree max(0, s — 3), 
with leading coefficient 

[k^-^]al{k-s) = til^ fors>3 (A.35) 

and first subleading coefficient 

[F-X(A;-s) = ^til^ fors>4 . (A.36) 

Proof. Using Lemma A. 3 we can rewrite (A. 34) as 

alik -s) = Ao,,,,.(3) + Ai,,,,,(3) , (A.37) 

where each of the two terms is, for integers k > s + 1, the restriction of a polynomial 
in k of degree s. Moreover, it follows from (A.20)/(A.21) that for a general value of 
A we have 

[F][Ao,.,fc(A) + Ai,,,fc(A)] = (A.38a) 

[k'-'][Ao,s,kW + A^,s,kW] = ^^~^iV (3-A) (A.38b) 

(s-1)! 

[F"2][Ao,.,fc(A) + A,,,fe(A)] = i0^4L^(3-A) (A.38c) 



2(s-2)! 

2 - \y-^ 

24(s-3)! 



[F-=^][Ao,.,fe(A) + Ai,,,fc(A)] = ^^1_!(3AV-15As2 + 18s2_a2s 

-19As + 66s - 24) (A.38d) 

[k'-'][Ao,s,kW + A^,,4\)] = -^l-^^(\h'-5\s' + 6s'-\h 

48(s — 4)! 

-19As + 66s - 24) (A.38e) 

Specializing to A = 3, we see that the first three coefficients vanish, while the next 
two coefficients are nonzero and given by (A.35)/(A.36). Therefore, af (/c — s) is the 
restriction to integers A; > s + 1 of a polynomial in k of degree exactly max(0, s — 3). 

D 
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A. 3 Proof of formula for a^{k — s) 

In this subsection we shall carry out an analogous analysis for the case of periodic 
transverse boundary conditions. It is convenient to prove first a technical lemma 
(analogous to Lemma A. 3) that deals with the outer sums we find in (3.63): 

Lemma A. 5 For any integers s > 1 and k > s, let us define 

B^AX) = (A; - s) V ^^ (^''~^~^] F.-p,s-i,p(A) , (A.39) 



where F^^a,b{X) is defined by (A. 7). We then have 

q\{s-q)\ 



Bs,kW = {k-s)J2 ^77^ Qi,s-iik - 9, A) , (A.40) 

g=0 "^ * 



where the polynomial Qe,a{b,X) is defined by (A. 9). In particular, Bsk{X) is the re- 
striction to integers k > s of a polynomial in k, A that is of degree s in k and A 
separately with leading coefficient 

[F]5,,fe(A) = -i(2-A)^ (A.41) 

and hence [k'^X'^] Bs^k{X) = {—l)''^~^^/s\. The next subleading coefficient is given by 

r-']B,4X) = ^(^^(2-Ar. (A.42) 

Proof. In (A.39) let us insert the result (A. 8) of Lemma A. 2 for -Fa;-p,s-i,p(A); we 
get 






-1)Mp-1)! 



p-l J {2p + s-k-l)\{k~p)\ 
xQk-p,s-iip,X) ■ (A.43) 

First of all, we see that for k = s, the sum is finite, and therefore Bk^k{X) = 0. So, 
we can restrict ourselves to the nontrivial cases A; > s + 1. Secondly, we notice, as 
in the proof of Lemma A. 3, that all terms with p < k — s — 1 vanish. So the only 
nonzero contributions to the sum (A.43) come from p > A; — s > L Using the variable 
q = k — p, we can rewrite (A.43) as 



R a^ (J. ^V-^g-^-A (-l)^(/c-g-l)! 

""-'^^^ = ^''-'^2^^[k-q-l)ik + S-2q-l)lql^--' 



ik-q,\) 



(I. ,^fk + s-2q-l\ (-1)'^+^ (A; - g - 1)! 

^"^-^^gi A;-g-l J (A; + .-2g-l)!g! ^^'-^^'-^'^) 



,=o^'(^-^)' 
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This expression for Bg^ki^) is manifestly a polynomial in k and A of degree at most s 
in each of the variables separately. Indeed, using (A. 10) and performing a binomial 
sum over q, we obtain (A. 41). 

The computation of the subleading term (A. 42) is similar to the one in Lemma A. 3. 
Using the expansion (A. 24) in the definition of Qq^s~i{k — q, A) that appears in (A. 44), 
we obtain the two leading terms: 



s-l 



Qq,s-i{k -g. A) 



E 



s-l 
m 



{k + s-2q-l 
\ rn, I 

m=0 

k^-^af + k'-\f + O {k'~') 



\S—l—m ^m \m 



g^A'"(A-r 



q—m 



where the coefficients a^ are given by 



.(2) 



(2) 

ay 



(x-irHs-i) 



s(A-l) 



+ g(2 - A) 



(A.45) 

(A.46a) 
(A.46b) 



The second task is to perform the sum over q in (A. 40). Because of the prefactor 
k—s, each coefficient a\ in (A. 46) gives a nonzero contribution to both [A;*~*"'"^]i?s^fc(A) 
and [k^~^]Bs,k{X). After some algebra we find 



[k']Bs,kW 



E 



'lY 



,(2)^ 



q\{s-q)\ 



q=0 

T.,o(l - A) 



(A.47a) 



r~']Bs,kW = Y. 



-lY 



5=0 

s(3 — s 



q\{s-q)\ 



.(2) 



sa 



(2) 



■T,,o(l - A) 



(^-1)(2-A) 
A-1 



T,,i(l-A) (A.47b) 



where the Tg^a are given in (A. 4). From these equations we can derive after some 
algebra the results (A.41)/(A.42). D 

Now we are ready to prove the main result of this subsection: 
Theorem A. 6 (= Proposition 3.7) For each integers > 1, the quantity 



T.'P 



a^ik-s) 



k 



p=l 



k — s fk 



k I ' ^^ ' p \ p — 1 



p — 1 



k—p 

E3' 

r=0 



k 



2p 



p 

k — p ~ r 

(A.48) 



is, when restricted to k > s, given by a polynomial in k of degree s, with leading 
coefficient 



[kTalik-s) 



■1) 



s+l 



(A.49) 
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and first subleading coefficient 

r-ralik-s) = ^^^^. (A.50) 

Furthermore, a^(0) = 0, so that the polynomial representing a^{k — s) for k > s has 
a factor k ~ s. 

Proof. The term 3'^( ^*) in (A. 48) vanishes for integers s > 1 and k > s, so we have 

alik-s) = B,,fe(3). (A.51) 

Therefore, the main assertions of Theorem A. 6 are an immediate consequence of 
Lemma A. 5. Finally, it is known from (3.65) that a|'(0) = for k > 1; and this 
corresponds to the case k = s, which lies within the regime k > s where a^{k — s) is 
given by a polynomial in k. D 



A. 4 Proof of formula for a^{k) 

Our last goal is to prove (3.43a). Let us start by proving two lemmas similar to 
Lemmas A. 3 and A. 5. 

Lemma A. 7 For any integer k > 1, let ai,fe(A) be defined as 

Then 

ai,.(A) = 1 + y]{-ir C'-^-P] A'^-^p-i . (A.53) 

U \ P J 

Proof. The inner sum in (A. 52), 

can be performed analogously to the proof of Lemma A. 2: If A; — 2p — 1 > 0, then 
( .^ ^) =0 for allr > k — 2p — I >0. Moreover, (^^^_ J = whenever k — p — r > 
p + 1 > 1, i.e. for all r < k — 2p — 1. Thus, when k — 2p — 1 > there is a single 
nonzero term coming from r = k — 2p — 1 and whose contribution is A'^'^^^^^. On the 
other hand, ii k — 2p — 1 < 0, then following similar steps to those of the proof of 
Lemma A. 2, we obtain a finite sum. Putting these two contributions together we get 

. ... ^ (p+1)! V Mr 

'"'^'^ ^ {2p-ky.f^^r\{k-p-ry.{2p + r + l-p) 

+A^-'^-'/[A;>2p+l], (A.55) 
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where I[A] is the indicator function of the event A, and in the first term it is not 
necessary to add I[k <2p + 1] because of the factor l/(2j9 — k)\. 
The outer sum in (A. 52) is easy: 



«..<^) - B-rf-ril^S;^ 



p=o X ^ / v-^ -^y ^^Q-{k-p-r)\{2p + r + l-p) 



k . 



^ ^ ^ ^ A'^-'P-i /[A; > 2p + 1] , (A.56) 



as the only nonzero contribution of the second term corresponds to p = k. (The 
contributions of the terms with < p < A; — 1 vanish because Euler's reflection 
formula.) We then obtain, 

p=0 \ P / 

Notice that the upper limit [_(k — 1)/2J comes from the condition k > 2p + 1; but we 
can replaced it by the simpler term /c — 1, as ( ^) = whenever p > k — l—p (i.e., 
2p > k — 1). We then obtain the claimed formula (A. 53). D 

Lemma A. 8 For any integer k > 1, let a2,fe(A) be defined as 

Then 

a2,kW = -1- (A.59) 

Proof. By inspection it is clear that a2,fc('^) = Ak{l,0;X) [cf. (A. 18)]. Then Lemma A. 3 
ensures that for any k > 1 

Afc(l,0;A) = -Qo,o(A;-l,A) = -1, (A.60) 

where we have also used (A. 9). D 

Now the strategy is to show that the quantities Ck{\) = ai^fc(A) + a2,fc(A) satisfy 
a second-order linear recursion, which reduces to that of the Fibonacci numbers F2k 
when A = 3. In particular we prove the following lemma: 

Lemma A. 9 For any integer k > 1, let Cfe(A) be defined as 

Cfc(A) = ai,fc(A)+a2,fc(A), (A.61) 

where the quantities aik{X) are defined in (A.52)/(A.58). Then, for any k > 3, the 
quantity 

/fe(A) = Cfc(A)-ACfe_i(A) + C,_2(A) (A.62) 

vanishes identically (i.e., /fc(A) = for all k > 3). 
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Proof. From Lemmas A. 7 and A. 8 it is clear that 

fc-i 



Ck{X) = E(-l)'(^ ^ ^)a^-^^-^ (A.63) 



p=0 

Then, the expression for /^(A) for any A; > 3 is given by 



P=0 \ P / p^Q \ P / 

+ J2i-ir\'-'^-'(^~^~^) . (A.64) 

p=0 \ P / 

The last term can be rewritten in a more convenient way if we change the summation 
variable to p + 1: 

This term is essentially the same as the second term in the sum (A.64), except for 
the binomial coefficients. Using the Pascal identity 

we obtain that 

P=0 \ P / p^Q \ P / 

= i-\)'-'5k,i = 0, (A.67) 

which is zero because we are assuming that k > 3. This completes the proof. D 

Now we are ready to prove the main result (3.43a): 
Theorem A. 10 For any integer k > 1 let us define the quantity [of. (3.34)] 

^i:MrC;!7)|-C?1(/-;-.)- <-) 

Then for any k > 1, a^{k) = F2k, where F2k is the {2k)-th Fibonacci number [30, 
sequences A001906/A088305]. 
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Proof. The quantity a^ik) can be written as 

al{k) = Cfc(3), (A.69) 

where Ck is given in Lemma A. 9 (see also Lemmas A. 7 and A. 8). For any A; > 3, the 
quantities a^(A;) satisfy the recurrence (A. 62) for A = 3, i.e. 

alik) - 3al_,{k - 1) + al_,{k - 2) = , (A.70) 

which is the same recurrence as the Fibonacci numbers F2n [30]. Furthermore, a 
direct computation shows that the initial conditions for the two sequences are the 
same: af (1) = 1 = F2 and 0^(2) = 3 = F4. Therefore, the sequences a^(A;) and F2k 
coincide. D 



B Upper zero-free interval for bipartite planar graphs 

Let G be a loopless planar graph. Then it is not hard to prove that PgW) > 
for all integers q > 5;^^ moreover, one of the most famous theorems of graph theory 
— the Four-Color Theorem [82-86] — asserts that Pg{<1) > holds in fact for all 
integers g > 4. 

It is natural to ask whether these results can be extended from integer q to real q. 
The answer is yes, at least in part: Birkhoff and Lewis [87] proved in 1946 that if G 
is a loopless planar graph, then Pg{(i) > for all real numbers q > 5.^^ Furthermore, 
they conjectured that Pg{.<i) > also for 4 < g < 5; and while no one has yet found 
a proof, no one has found a counterexample either, so it seems plausible (in the light 
of the Four-Color Theorem) that the conjecture is true. 

Now, some planar graphs can be colored with three or even two colors; this means 
that their chromatic polynomials Pcio) are strictly positive for integers g > 3 or 
g > 2, respectively. Can these bounds can be extended to real g? That is, if G is a k- 
colorable planar graph, do we have Pg{(i) > for all real q>kl Woodall [88, p. 142] 
conjectured that the answer is yes. For /c = 4, this is the conjecture of Birkhoff 
and Lewis mentioned above. For /c = 3, however, Thomassen [89, pp. 505-506] 
has shown that Woodall's conjecture is false: there exist 3-colorable (and in fact 
2-degenerate^^) planar graphs with real chromatic roots greater than 3. Indeed, by 



■^"'This is the Five-Color Theorem, which goes back to Heawood in 1890. For a proof, see e.g. [81, 
Theorem V.8, pp. 154-155]; or for an elegant alternate proof of an even stronger result, see [81, 
Theorem V.12, pp. 161-163]. 

^^See also Woodall [88, Theorem 1] and Thomassen [89, Theorem 3.1 ff.] for alternate proofs 
of a more general result. These theorems are, in turn, consequences of a stronger and even more 
general result for matroids that was proved but not stated (!) twenty years earlier by Oxley [90]; see 
Jackson [91, Theorem 38]. 

36^ graph G is said to be k-degenerate if every subgraph H C G has at least one vertex of 
degree < k. It is not difhcult to prove [81, Theorem V.l, p. 148] that every fc-degenerate graph is 
(fc + l)-colorable. 
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combining Thomassen's construction [89, proof of Theorem 3.9] with Royle's [92] 
recent construction of planar graphs with real chromatic roots arbitrarily close to 4, 
we see that there exist 3-colorable (and in fact 2-degenerate) planar graphs with real 
chromatic roots arbitrarily close to 4. 

In Ref. [6] we showed that Woodall's conjecture is false also for k = 2: there 
exist 2-colorable (i.e. bipartite) planar graphs with real chromatic roots greater than 
2. For example, the 4p x 6f square lattice has chromatic roots at g ~ 2.009978 and 
q ^ 2.168344. For the cases 8f x np and 8p x np, we also observed numerically [6] that 
there are real chromatic roots tending to B^ = {3 + V5)/2 ^ 2.618034 from below as 
n ^ oo. This led us to modify Woodall's conjecture as follows [6, Conjecture 7.5]: 

Conjecture B.l [6] Let G be a bipartite planar graph. Then Pg{q) > for real 
q>B, = i3 + V5)/2. 

However, the numerical results reported in the present paper (see Tables 15 and 
16) now show that this conjecture is also false. Indeed, it appears that all strip graphs 
9f X np (with odd n > 39) and lip x np (with odd n > 23) have two real roots greater 
than B^. 

So, not only is Conjecture B.l false; it is actually false for two distinct reasons. 
First of all, the real roots converging to the isolated limiting point B^ need not do so 
only from below. (Empirically we find that, for m = 9 and 11, these roots converge to 
i?5 with parity (—1)"+^.) Secondly, real roots can converge to the real crossing point 
qo{m) > Br, when one exists. ^"^ (Empirically we find that such a crossing point exists 
at least for odd widths m and that the real roots converging to it exist for sufficiently 
large odd lengths n.) 

We expect both of these behaviors to persist for all larger odd widths m. In 
particular, since we expect that qo{m) f 3 as m ^ oo (see Section 6.2), we expect that 
m-p X np square-lattice strips with m, n both odd and large can have real chromatic 
roots arbitrarily close to 3. However, in most cases of which we are aware [6-10,92], 
the convergence of real chromatic roots to qo{fn) occurs only from below. ^^ Therefore, 
the following weakened version of Conjecture B.l is plausible: 



■^^A similar phenomenon was exploited recently by Royle [92] to provide examples of plane trian- 
gulations with real chromatic roots converging to 4 from below. His graphs are 4p x np triangular 
lattices with carefully chosen endgraphs adjoined at top and bottom (sec [93] for a general analysis 
of such strips- with-endgraphs). The key fact underlying this construction is that go(tri,4p) ~ 4. 

^^The only exceptions we know are as follows: 

(a) For square-lattice strips of width m = 2 with cyclic boundary conditions [94] , and triangular- 
lattice strips of width m = 2,3 with cyclic and toroidal boundary conditions [9, 10], it appears 
that there exists an infinite sequence of lengths n with real zeros converging to qo{m) from 
above. 

(b) For triangular-lattice strips of width m = 4, 5 with cyclic boundary conditions [9] and width 
m ~ 5,7 with toroidal boundary conditions [10], and for square-lattice strips of width m = 3 
with cyclic boundary conditions [9], there are real zeros q > qo{m) for some small lengths n, 
but the ultimate convergence of real zeros to qoi'm) appears to be only from below. 

There are also some cases in which the maximal real zero is exactly equal to qo {m) . Finally, in some 
cases (e.g., triangular-lattice strips of width m = 4,5,6,7 with toroidal boundary conditions [10]) 
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Conjecture B.2 Let G be a bipartite planar graph. Then Pcio) > for real q > 3. 

(Note that the honeycomb lattice, which is bipartite and planar, has Qc = B^ < 3.) 

We can also give some plausible "physics" considerations that provide additional 
support for a slightly weakened version of Conjecture B.2. Some two-dimensional 
antiferromagnetic models at zero temperature have the remarkable property that they 
can be mapped onto a "height" (or "interface" or "SOS-type") model (see e.g. [74] and 
the references cited there). Experience tells us that when such a representation exists, 
the corresponding zero-temperature spin model is most often critical. ^^ In particular, 
when the g-state zero-temperature Potts antiferromagnet on a lattice C admits a 
height representation, one expects that q = qc{C)^^ This prediction is confirmed in all 
heretofore-studied cases: 3-state square-lattice [71,73,74,96], 3-state kagome [97,98], 
4-state triangular [99], and 4-state on the line graph (= covering lattice) of the square 
lattice [98,100]. Now, the height representation of the square-lattice zero-temperature 
3-state Potts antiferromagnet (i.e., the square-lattice chromatic polynomial at g = 3), 
as presented in [74], generalizes immediately from the square lattice to an arbitrary 
finite or infinite plane quadrangulation G. This suggests that qdjC) = 3 for all infinite 
periodic lattices that are plane quadrangulations. It in turn provides some support 
for the following weakened version of Conjecture B.2: 

Conjecture B.3 Let G be a plane quadrangulation. Then Paiq) > for real q > 3. 

Surprisingly, however, we know of a potential counterexample to Conjectures B.2 
and B.3! Consider the diced lattice [101, Figure 15], which is the dual of the kagome 
lattice and is a plane quadrangulation. Even though there exists a height representa- 
tion for the 3-state Potts antiferromagnet on the diced lattice, this model is in a non- 
critical ordered phase at zero temperature. Indeed, Feldmann et al. [102] found (by 
dualizing the numerical data of Ref. [103] for the kagome-lattice Potts model) that the 
3-state diced-lattice Potts model has a critical point at Vc = —0.8607 ±0.0008, which 
lies at nonzero temperature in the physical antiferromagnetic region — 1 < w < 0. 
(By contrast, the 4-state model was found to have its antiferromagnetic critical point 
in the unphysical region at Vc = — 1.18 ± 0.02, so that the 4-state antiferromagnet 



there are sequences of complex zeros converging to qo{m) but apparently no sequences of real zeros 
converging to qo{'m). 

Let us remark that the tables of real chromatic roots in Refs. [6-9] show only the cases of strips 
m X (km) with k ~ 1, . . . , 10. Moreover, in Ref. [10] such tables are absent. Therefore, we have here 
used the transfer matrices reported in those papers to recompute all the chromatic polynomials for 
strips m X n with 1 < n < lOrn (and in some cases for even larger lengths) and their corresponding 
chromatic roots. 

■^''Some exceptions are the constrained square-lattice 4-statc antiferromagnetic Potts model [73] 
and the triangular-lattice antiferromagnetic spin-s Ising model for large enough s [95] , both of which 
appear to lie in a non-critical ordered phase at zero temperature. 

^°Here qdC) is the value (which is conjectured to exist) such that for q > qdC) the antiferro- 
magnetic Potts model has exponential decay of correlations uniformly at all temperatures, including 
zero temperature, while for q — qdC) the model has a zcro-tcmpcrature critical point. 



is disordered for all temperatures, including zero temperature.) We have very re- 
cently [79] confirmed this scenario by Monte Carlo simulations: for g = 3 we find a 
critical point at Vc = —0.860599 ±0.000004, while for g = 4 we find a finite correlation 
length uniformly in the region —1 < f < 0. Furthermore, for g = 3 we are able to 
give a mathematically rigorous proof of the existence of a phase transition at nonzero 
temperature, using a (computer-assisted) Peierls argument. These results indicate 
that 3 < gc(diced) < 4, contrary to our prediction that gc = 3 for all (regular) plane 
quadrangulations. Indeed, linear interpolation from the predictions of [102] would 
suggest gc( diced) ~ 3.44. 

Now, one strongly expects that finite pieces of the diced lattice with free or 
cylindrical boundary conditions will have chromatic roots tending to gc (diced) in 
the infinite-volume limit. What is less clear, however, is whether the roots converging 
to Qc will necessarily be real roots. If they are, then a sufficiently large finite piece 
of the diced lattice, with boundary conditions arranged so that the graph is a plane 
quadr angulation (or at least planar and bipartite), would provide a counterexample 
to Conjecture B.3 (or at least to Conjecture B.2). But if they are not, then the 
fact that gc > 3 is in no way incompatible with Conjectures B.2 and B.3. Please 
note also that the roots converging to Qc can be complex for some choices of end- 
graphs and real for others; this possibility was recently exploited by Royle [92] to 
provide examples of plane triangulations (based on 4p x np triangular-lattice strips 
with suitable endgraphs [93]) having real chromatic roots converging from below to 
gc(tri) = 4. It would be very interesting to undertake a transfer-matrix study of the 
diced-lattice chromatic polynomial along the lines of [6-10,92,93], in order to test 
whether the roots converging to Qc can be real, at least for some choices of planar 
bipartite endgraphs. 
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Table 1: Coefficients bl{m) of tlie large-g expansion of tlie dominant eigenvalue Af for free boundary conditions. For each 
1 < -m < 12, we include all coefficients &f (m) up to A; = 15. For the whole data set up to A; = 40, see the Mathematica file 
data_FREE.m included in the on-line version of the paper at arXiv.org. Those data points below the stair-case-like line satisfy 
^ > <i„(^) [cf. (4.11)]. 
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Table 2: 



Coefficients k4{'m') 



of the large-g expansion of log(g ""A^) where Af is the dominant eigenvalue for free boundary 



conditions. For each 1 < m < 12, we include all coefficients ^i"^) up to A; = 15. For the whole data set up to A; = 40, see 
the Mathematica ffie data_FREE.m included in the on-line version of the paper at arXiv.org. Those data points below the 
stair-case-like line satisfy m > m^^^i^k) [cf. (4.11)]. 
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Table 3: Coefficients A^{m) for free boundary conditions [cf. (4.16)] for 6 < A; < 21 and 1 < -m < 12. For < A; < 6 and 
1 < m < 12, the coefficients A|(m) vanish. For the whole data set up to /c = 33, see the Mathematica file data_FREE.m 
included in the on-line version of the paper at arXiv.org. Those data points below the stair- case- like line satisfy m > Tn^^^^{k) 
[cf. (4.11)] and are therefore zero. 
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Table 4: Coefficients h\{m) of the large-g expansion of the dominant eigenvalue A^ for cylindrical boundary conditions. For 
each 1 < -m < 13, we include all coefficients h\{m) up to A; = 15. For the whole data set up to A; = 40, see the Mathematica 
ffie data_CYL.m included in the on-line version of the paper at arXiv.org. Those data points below the stair- case- like line satisfy 
m > 'n^min(^) = k + 2. 
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Table 5: Coefficients A;c^(m) of the large-g expansion of log(g~"^A^) where A^ is the dominant eigenvalue for cylindrical boundary 
conditions. For each 2 < -m < 13, we include all coefficients c^(m) up to A; = 15. For the whole data set up to /c = 40, see 
the Mathematica file data_CYL.m included in the on-line version of the paper at arXiv.org. Those data points below the 
stair-case-like line satisfy m > fn^^^^lk) = k + 2. 
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Table 6: Coefficients b^{m) of the large-g expansion of the dominant eigenvalue A^ for cylindrical boundary conditions. For 
each 1 < -m < 13, we include all coefficients &f (m) up to A; = 15. For the whole data set up to A; = 40, see the Mathematica 
file data_CYL.m included in the on-line version of the paper at arXiv.org. Those data points below the stair-case-like line satisfy 
m > k. 
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Coefficients kcfi 



defined in (4.10). For each 2 < m < 13, we include all coefficients c^{m) up to k = 15. For the 



Table 7: uoetticients fcc^ [m^ 

whole data set up to A; = 40, see the Mathematica file data_CYL.m included in the on-line version of the paper at arXiv.org. 

Those data points below the stair-case-like line satisfy m > 'fn^^^^{k) = max(A;, 2). 
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Table 8: Modified coefficients A^{m) for cylindrical boundary conditions for 2 < k < 16 and 2 < m < 13. For k = 1 and 
2 < -m < 13, the coefficients A^{m) vanish. Those data points below the stair-case- like line satisfy m > 'fn^^:^j^{k) = niax(A;, 2) 
and are therefore zero. 
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12262136381593 


-38841312202313 


45 


5101697398582 


-24192583755347 


104412348649015 


46 


12146149238921 


-118764516172484 


448320847685638 


47 


-598931311074 


-130312353695974 





Table 9: Large-g series expansions for the bulk, surface and corner free energies of 
the square-lattice zero-temperature Potts antiferro magnet, in terms of the variable 



l/(g — 1). The coefficients are defined by e-'^''""' = [(g — l)^/q] Xlfclo ^kZ^, 



Djsuri 



Yll^=oPkZ^ and 6-^'=°™'="' = YlT={ilkz''- The coefficients a^ for /c < 36 were obtained in 
Ref. [12]; all the other coefficients are new. The coefficient fi^-j is conjectural. 
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CO 



A^ 


Z\ 


^2 


^3 


30 


-0.27(1) ± 0.773(7) i 
0.6(2) 




0.502(3) ± 0.292(5) i 
0.2(2) 


40 


-0.2811(5) ± 0.7757(7) i 
0.47(3) 


0.45(2) ± 0.53(1) i 
0.4(5) 


0.501(2) ± 0.289(2) i 
0.4(1) 


50 


-0.2811(2) ± 0.7752(1) i 
0.50(1) 


0.452(6) ± 0.538(5) i 
0.3(2) 


0.49999(2) ±0.28868(2) 2 
0.501(5) 


60 


-0.28112(3) ± 0.77520(2) i 
0.500(3) 


0.4480(9) ± 0.538(1) i 
0.48(7) 


0.499997(5) ± 0.288676(6) i 
0.501(2) 


70 


-0.281118(3) ± 0.775201(4) i 
0.4999(6) 


0.44778(4) ± 0.53827(8) i 
0.500(4) 


0.499999(2) ± 0.288675(1) i 
0.5001(5) 


80 


-0.281118(2) ± 0.775201(2) i 
0.5000(3) 


0.447783(5) ± 0.538251(8) i 
0.5000(7) 


0.500000(1) ± 0.2886751(1) i 
0.5001(4) 


90 


-0.281118(2) ± 0.775201(2) i 
0.5000(3) 


0.447784(1) ± 0.538250(2) i 
0.4999(3) 


0.4999997(9) ± 0.2886752(8) i 
0.5001(3) 


100 


-0.2811174(5) ± 0.7752009(6) i 
0.5000(1) 


0.447784(1) ± 0.538249(2) i 
0.5000(2) 


0.4999998(8) ± 0.2886752(8) i 
0.5000(2) 


Exact 


-0.2811172691 ± 0.7752009092 i 
0.5 


0.44778393657 ± 0.5382490441 i 
0.5 


0.5 ± 0.2886751346 i 
0.5 


^i 


0.8245989138 


0.7001589015 


0.5773502692 



Table 10: Estimates of the position and exponent of the "correct" singularities for the square- lattice strip of width L = 3 and 
free boundary conditions. The column A^ indicates the maximum order of the series in z = l/(g — 1) that was used in the 
analysis. The row marked 'Exact' shows the known exact results. The last row (marked \zi\) shows the absolute value of the 
corresponding exact singularity. Blank entries mean that we were unable to obtain the corresponding estimates. 



N 


Z4 


Z5 


Z6 


Z7 


30 


0.36(1) ± 0.66(2) i 
0.4(5) 


-1.000(5) 
-0.1(2) 






40 


0.368(7) ± 0.654(8) i 
0.0(2) 


-1.0000(1) 
-0.001(6) 






50 


0.3695(6) ± 0.6524(1) i 
-0.03(5) 


-1.0000(1) 
0.000(6) 






60 


0.3695(2) ± 0.6516(4) i 
-0.01(3) 


-1.00000(6) 
0.000(4) 






70 


0.3695(1) ± 0.6514(2) i 
0.00(1) 


-1.00000(2) 
0.000(2) 




0.487(5) ±0.427(4)2 
2.03(9) 


80 


0.36947(2) ± 0.65136(2) i 
-0.001(3) 


-1.00000(2) 
0.000(2) 


0.59(5) 
0.0(6) 


0.488(2) ± 0.428(1) i 
2.01(4) 


90 


0.369460(3) ± 0.651363(2) i 
0.0001(6) 


-1.00000(2) 
0.000(2) 


0.595(1) 
0.001(1) 


0.4878(5) ± 0.4277(2) i 
2.000(2) 


100 


0.369460(3) ± 0.651363(3) i 
0.0000(6) 


-1.000000(9) 
0.0000(9) 


0.5945(4) 
0.000(8) 


0.4877(4) ± 0.4277(2) i 
2.000(2) 



Table 11: Estimates of the position and exponent of the spurious singularities for 
the square-lattice strip of width L = 3 and free boundary conditions. The column 
N indicates the maximum order of the series in z = l/{q — 1) that was used in 
the analysis. Blank entries mean that we were unable to obtain the corresponding 
estimates. 
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CO 
CO 



N 


21 


22 


23 


24 


25 


30 


0.331(2) ± 0.238(2) J 
0.5(2) 


-0.068(1) ± 0.4795(8) i 
0.53(8) 


0.29(2) ± 0.49(3) J 
-0.3(3) 






40 


0.3315(4) ± 0.2370(4) i 
0.53(8) 


-0.0682(2) ± 0.4798(2) i 
0.51(3) 


0.30(2) ± 0.48(2) i 
-0.4(6) 






50 


0.33162(6) ± 0.23712(5) i 
0.50(2) 


-0.068171(3) ± 0.479859(4) i 
0.5005(9) 


0.30(1) ± 0.48(3) J 
-0.2(6) 






60 


0.331635(2) ±0.237117(2)1 
0.4995(9) 


-0.0681711(4) ± 0.4798607(9) i 
0.5001(3) 


0.30(2) ± 0.46(3) J 
0.3(13) 






70 


0.3316353(8) ± 0.237116(1) i 
0.4998(4) 


-0.0681712(2) ± 0.4798609(5) i 
0.5000(1) 


0.29(2) ± 0.46(3) J 
0.4(12) 






80 


0.3316354(7) ± 0.2371156(1) i 
0.4998(4) 


-0.0681712(2) ± 0.4798609(2) i 
0.50000(6) 


0.29(1) ± 0.45(2) J 
0.6(10) 


0.76(3) 
-0.9(6) 




90 


0.3316354(6) ± 0.2371155(9) i 
0.4999(3) 


-0.0681712(2) ± 0.4798609(2) i 
0.50001(7) 


0.28144(4) ± 0.45208(5) i 
0.503(6) 


0.75(1) 
-0.7(4) 




100 


0.3316354(7) ± 0.237116(1) i 
0.4999(4) 


-0.0681712(1) ± 0.4798609(3) i 
0.50000(9) 


0.281427(4) ±0.452108(3)1 
0.5000(4) 


0.754(9) 
-0.7(3) 




Exact 


0.3316354418 ± 0.2371152471 i 
0.5 


-0.0681712693 ± 0.4798609413 i 
0.5 


0.2814289723 ± 0.4521062477 i 
0.5 


0.7398155434 
0.5 


0.7978491474 
0.5 


kd 


0.4076833412 


0.4846791154 


0.5325432618 


0.7398155434 


0.7978491474 



Table 12: Estimates of the position and exponent of the "correct" singularities for the square-lattice strip of width L = 4 and 
periodic boundary conditions. The column A^ indicates the maximum order of the series in z = l/(g — 1) that was used in the 
analysis. The row marked 'Exact' shows the known exact results. The last row (marked \zi\) shows the absolute value of the 
corresponding exact singularity. 



N 


Z6 


Z7 


30 


-0.999(4) 
0.0(1) 




40 


-1.000(1) 
0.00(3) 




50 


-1.00000(4) 
0.000(2) 


0.30(2) ± 0.47(2) i 
0.3(12) 


60 


-1.00000(1) 
0.0002(9) 


0.29(1) ± 0.46(2) i 
0.5(10) 


70 


-0.999999(9) 
0.0001(5) 


0.29(2) ± 0.46(3) i 
0.3(10) 


80 


-1.00000(1) 
0.0001(7) 


0.286(2) ± 0.517(2) i 
0.0(1) 


90 


-1.000000(2) 
0.0000(2) 


0.2854(4) ± 0.5180(4) i 
0.00(3) 


100 


-1.000000(4) 
0.0000(4) 


0.28528(4) ± 0.51814(6) i 
0.000(4) 



Table 13: Estimates of the position and exponent of the spurious singularities for the 
square-lattice strip of width L = 4 and periodic boundary conditions. The column 
N indicates the maximum order of the series in z = l/(g — 1) that was used in 
the analysis. Blank entries mean that we were unable to obtain the corresponding 
estimates. 



100 



mp 


Re go 


Imgo 


Type 


3f 


2 





S 


4f 


2.2283590792 





+ 


5f 


2.4284379020 





s 


6f 


2.5286467909 





+ 


7f 


2.6062482130 





s 


8f 


2.6602596967 


0.0007413717 


D 


9f 


2.7016599568 





s 


IOf 


2.7343903604 


0.0003924978 


D 


11f 


2.7608973951 





S 


12f 


2.782817590 


0.00018700 


D 



Table 14: Values of qo{iTi) for square-lattice strips with free boundary conditions. 
Type "S" means that the limiting curve Bm crosses the real axis at a single point, 
"D" means that the limiting curve does not cross the real g-axis but has two complex- 
conjugate endpoints nearby, and "+" means that the limiting curve contains a seg- 
ment on the real g-axis (we define go to be the lower end of that segment). For points 
of type "D", we also show the imaginary part. The data for m < 7 are taken from [6]. 
Notice that the value for m = 8 reported in [6] is wrong; the correct value is displayed 
here. 
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n 


4th Zero - B^ 


5th Zero 


30 


-8.9459389900 x lO^^ 




32 


-5.6487223660 x lO^^ 




34 


-3.1962919070 x 10-^ 




36 


-1.5712164060 x lO'^ 




38 


-6.6617400470 x 10"^ 




39 


5.8575529850 x 10"^ 


2.6298556468 


40 


-2.5161135990 x 10"^ 




41 


1.7270275260 x 10"^ 


2.6345567384 


42 


-8.8890257250 x 10"^ 




43 


5.5346483490 x 10"^ 


2.6384105545 


44 


-3.0434570840 x 10"^ 




45 


1.8199961300 x 10"^ 


2.6417380816 


46 


-1.0281765140 x 10"^ 




47 


6.0420071920 x 10"^ 


2.6446771074 


48 


-3.4539716110 X 10"^ 




49 


2.0135551740 x 10"^ 


2.6473078802 


50 


-1.1574596530 x 10"^ 




51 


6.7216127470 x 10"^ 


2.6496845069 


52 


-3.8744221380 x 10"^ 




53 


2.2455458550 x 10"^ 


2.6518466022 


54 


-1.2962041990 x 10"^ 




55 


7.5047689080 x 10"^ 


2.6538246782 


56 


-4.3353170180 x 10"^ 




57 


2.5086443380 x 10"^ 


2.6556430336 


58 


-1.4497910190 x 10"^ 




59 


8.3866240660 x lO'^ 


2.6573214747 


60 


-4.8479232330 x 10"^ 




61 


2.8038883660 x 10-^ 


2.6588764247 


62 


-1.6210151380 X 10-9 




63 


9.3745081850 x lO^^" 


2.6603216811 


64 


-5.4201053100 x 10-^° 




65 


3.1343270220 x 10-^° 


2.6616689582 


66 


-1.8122675760 x 10-^° 




67 


1.0479599710 x 10-^° 


2.6629282884 


68 


-6.0594535690 x lO-^^ 




69 


3.5038676450 x lO^i^ 


2.6641083276 


70 


-2.0260148680 x lO^^^ 




oo 





2.7016599568 



Table 15: Fourth and fifth real chromatic zeros for a strip 9f x n-p with free boundary 
conditions. As the fourth zero converges rapidly to Br,, we show their difference. The 
last row (labelled with n = oo) shows the infinite-length limit; in particular, the value 
ofgo(9). 
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n 


4th Zero - B^ 


5th Zero 


14 


-3.6516138320 x 10"^ 




16 


-1.9375298640 x 10"^ 




18 


-8.1115634110 X 10-3 




20 


-2.2141892360 x lO^^ 




22 


-3.6952234840 x 10"^ 




23 


2.7283542170 x 10"^ 


2.6360127941 


24 


-4.7440730960 x 10"^ 




25 


2.3461500250 x 10"^ 


2.6467690175 


26 


-5.5660625770 x 10"^ 




27 


2.2973640450 x 10^^ 


2.6554417441 


28 


-6.2842068380 x 10"^ 




29 


2.3500000310 x 10"^ 


2.6628475620 


30 


-6.9454825890 x 10^^ 




31 


2.4549166270 x 10'^ 


2.6693333011 


32 


-7.5798800910 x 10"^ 




33 


2.5935418370 x 10"^ 


2.6750988717 


34 


-8.2103756740 x 10-^° 




35 


2.7571779370 x lO-^o 


2.6802756665 


36 


-8.8542077070 x lO-^^ 




37 


2.9415010300 x lO^^^ 


2.6849567990 


38 


-9.5240493620 x lO^^^ 




39 


3.1444559750 x lO^^^ 


2.6892116307 


40 


-1.0229343110 x lO^^^ 




41 


3.3652761350 x lO^^^ 


2.6930937453 


42 


-1.0977435000 x lO'^^ 




43 


3.6039742950 x lO'^^ 


2.6966457592 


44 


-1.1774405610 x lO^^^ 




45 


3.8610612020 x lO^^^ 


2.6999024444 


46 


-1.2625636940 x lO'^^ 




47 


4.1373846520 x 10"^^ 


2.7028928687 


48 


-1.3536187010 x lO'^^ 




49 


4.4340363840 x lO^^^ 


2.7056419158 


50 


-1.4511035420 x lO^^^ 




51 


4.7522986580 x lO'^^ 


2.7081713766 


52 


-1.5555245880 x 10^^^ 




53 


5.0936147060 x lO'^^ 


2.7105007206 


54 


-1.6674076830 x lO^^^ 




oo 





2.7608973951 



Table 16: Fourth and fifth real chromatic zeros for a strip lip x np with free boundary 
conditions. As the third zero converges rapidly to B^, we show their difference. The 
last row (labelled with n = oo) shows the infinite-length limit; in particular, the value 
ofgo(ll). 
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DA method K=2 for fg^: N=50 



in 
d 

I 



— 1 — I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — r— 



in 
6 



O i-Zs 




J I I I I I I I I I I I I I I I I L 



-1 



-0.5 

Re z 

fa) 
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DA method K=2 for fgp: N=100 



d 

I 



— 1 — I — I — I — I — I — I — I — I — I — I — I — I — r—n — I — I — r- 



d 



o i-zs 




_i I I I I I I I I I I I I I n I I i_ 



-1 



-0.5 




Re z 

(b) 
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Figure 1: Histograms of the non-defective zeros for the square lattice strip of width 
L = 3p using K = 2 differential approximants with A^ = 50 (a) and N = 100 (b) 
coefficients of the free-energy series expansion. Each cell is a square of size 0.04, and 
the bin count is indicated with a gray-scale code: there are 16 gray tones ranging 
from black (largest counts) to white (smallest counts). We also show the limiting 
curve B^p. 
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Sq lattice L=4p 



■1 -0.5 0.5 

Re z = Re l/(q-l) 



w o 



a o - 






-"I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — I — r 




■1 -0.5 0.5 

Re z = Re l/(q-l) 

(b) 



Figure 2: Plot of the limiting curves for the square-lattice strips of widths 3f (a) and 
4p (b) in the complex z = l/{q ~1) plane. We also show the location of the physical 
singularities (i.e., the endpoints of the limiting curve, depicted in blue D) and of the 
spurious ones (red o). The dashed line shows the circle with the radius of convergence 



obtained by a raw fit using Tc 
to those of the tables. 



lim inf ^ 



-i/n rjnj^g point labels correspond 
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DA method K=l for exp(f, ,^) 



DA method K=2 for exp(f, ,^) 



id 
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id 
6 



LO 
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Figure 3: Histograms of the non-defective zeros for the bulk-free-energy series (5.23) 
using K = 1 {a) and K = 2 (h) differential approximants with N = 47 coefficients. 
Gray-scale codes are as in Figure 1. For comparison, we show the limiting curve Bm 
for the square-lattice strip of width m = 7 with toroidal boundary conditions [10]. 
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DA method K=2 for f 



bulk 



r 



o 



N 



o 

I 



T — I — r 



Ik. 



T 



T — I — r 




J I L 



1 



T 
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Re z 
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Figure 4: Histogram of the non-defective zeros for the unexponentiated bulk-free- 
energy series (5.25) using K = 2 differential approximants with N = Al coefficients. 
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DA method K=l for exp(f ,) 



DA method K=2 for exp(f ,) 



d 



6 

I 




Figure 5: Histograms of the non-defective zeros for the surface-free-energy series (5.33) 
using i^ = 1 (a) and K = 2 (h) differential approximants with N = 47 coefficients. 
Codes are as in Figure 1. We show the hmiting curve Bm for the square-lattice strip 
of width L = 11 with toroidal boundary conditions as a solid black curve. 



108 



DA method K=l for exp(f ) 

^ ^ corner' 



DA method K=2 for exp(f ) 

^ ^ corner' 
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DA method K=2 for f. 



— I — I — I — r- 
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Figure 6: Histograms of the non-defective zeros for the corner-free-energy series (5.32) 
using K = 1 (a) and K = 2 (h) differential approximants with A^ = 46 coefficients. 
Panel (c) shows the analogous K = 2 histogram for the unexponentiated series (5.34). 
Codes are as in Figure 1. We show the limiting curve Bm for the square-lattice strip 
of width L = 11 with toroidal boundary conditions as a solid black curve. 
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Bulk free energy 



Surface free energy 




■1 -0.5 0.5 

Re z = Re l/(q-l) 




-1 -0.5 0.5 

Re z = Re l/(q-l) 

(b) 



Corner free energy 




-0.5 
Re z 



0.5 

Re l/(q-l) 



Figure 7: Singularities for the bulk (a), surface (b), and corner (c( free energies in 
the complex plane oi z = l/(g — 1). In each panel, we show the limiting curve 
Bm for the square-lattice strip of width L = 7 (resp. L = 11) with toroidal (resp. 
cylindrical) boundary conditions as a solid black (resp. dashed green) curve. We 
also show the location of the singularities found in the series analysis (blue D). The 
dashed black line shows the circle with the radius of convergence estimated from 



liminf„^„o \an 



-l/n 
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Zeros sq lattice L= 9^ 



Zeros sq lattice L= lOp 



CV2 - 
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Figure 8: Zeros of the chromatic polynomial for the square lattices of widths 9f (a), 
lOp (b), lip (c), and 12f (d). For each width L, we show the chromatic zeros of the 
strips Lp X (5L)f (□ black) and Lp x (10L)f (o red). 
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